2, 2008 17:43 WSPC/INSTRUCTION FILE Yang-Lee 



International Journal of Modern Physics B 
© World Scientific Publishing Company 



STATISTICAL MECHANICS OF EQUILIBRIUM AND 
NONEQUILIBRIUM PHASE TRANSITIONS: 
THE YANG-LEE FORMALISM 



lOANA BENA 
Theoretical Physics Department, University of Geneva, 
24, Quai Ernest Ansermet, GH-1211 Geneva 4, Switzerland 
loana. Bena@physics.unige. ch 

MICHEL DROZ 

Theoretical Physics Department, University of Geneva, 
24, Quai Ernest Ansermet, GH-1211 Geneva 4, Switzerland 
Michel. Droz&physics.unige. ch 

ADAM LIPOWSKI 

Quantum Physics Division, Faculty of Physics, A. Mickiewicz University 
Ul. Umultowska 85, 6I-6I4 Poznan, Poland 
lipowski@amu. edu.pl 

Received Day Month Year 
Revised Day Month Year 

Showing that the location of the zeros of the partition function can be used to study 
phase transitions, Yang and Lee initiated an ambitious and very fruitful approach. We 
give an overview of the results obtained using this approach. After an elementary intro- 
duction to the Yang-Lee formalism, we summarize results concerning equilibrium phase 
transitions. We also describe recent attempts and breakthroughs in extending this theory 
to nonequilibrium phase transitions. 
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1. Introduction 

Boiling of water, denaturation of a protein, or formation of a percolation cluster in a 
random graph are examples of phase transitions. Description and understanding of 
these ubiquitous phenomena, that appear in physics, biology, or chemistry, remains 
one of the major challenges of thermodynamics and statistical mechanics. 

The simplest and relatively well-understood phase transitions appear in sys- 
tems in thermal equilibrium with their environment. A prototypical example of 
such an equilibrium phase transition is the liquid-gas transition. Crossing the coex- 
istence curve, by varying e.g., the temperature, one can move between vapor and 
liquid phases. Such a transition is called first-order (or discontinuous) because in 
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the framework of thermodynamics this phenomenon is understood in terms of ther- 
modynamic potentials showing discontinuities in their first derivatives with respect 

to the control parameters. The coexistence line terminates at a particular point, 
the critical point, where both phases become undistinguishable. Close to this point, 
second derivatives of the thermodynamic potential exhibit power-law singularities 
and this indicates a second-order (or continuous) phase transition. In the late fifties, 
giving a microscopic description of first- and second-order phase transitions became 
a challenge, that culminated in the development of scaling theories and then renor- 
malization group methods. Although these techniques are computationally very 
effective and reliable, they do not explain all the aspects of phase transitions. In 
particular, it is not clear how free energy can develop singularities at the first-order 
phase transitions. To address this problem, Yang and Lee developed an approach 
where the partition function (that is used to calculate the thermodynamic poten- 
tials) is considered as a function of complex control parameters. Singularities of 
the thermodynamic potentials, given as zeros of the partition function, were then 
shown to accumulate exactly at the transition point. This approach, that was later 
on generalized and extended to various systems, provides a lot of information about 
equilibrium phase transitions. 

However, an equilibrium system is rather often only an abstract idealization. 
Typically, a given system is not in equilibrium with its environment, but is ex- 
changing matter and/or energy with it; fluxes are present in the system. The usual 
way to describe the physics of such a system is to write down a master equation 
for the time-dependent probability that a given mesoscopic state is realized at time 
t. The physics is then embedded in the transition rates between these mesoscopic 
states. In the long-time limit the system eventually settles into a steady-state and, 
depending on the values of some control parameters or the initial configuration, 
different nonequilibrium steady states (NESS) can be reached, and sometimes a 
breaking of ergodicity takes place. Upon varying some control parameters the sys- 
tem may change its steady state and this is called a nonequilibrium phase transition. 
In the particular situation in which the so-called detailed balance condition is ful- 
filled, one recovers the case of equilibrium phase transitions. Our understanding 
of NESS from a microscopic point of view is not as advanced as for equilibrium. 
Accordingly, it is legitimate to try to extend the Yang-Lee theory to NESS. It was 
argued recently, and illustrated on several examples, that Yang-Lee ideas should 
apply, at least to some extent, to nonequilibrium phase transitions. 

This review provides an overview of the Yang-Lee theory for equilibrium phase 
transition and describe recent advances in applying it to nonequilibrium systems. 
In addition, we included an extensive bibliographical list that offers an interested 
reader the relevant references where one may find further details on a more specific 
subject. Our presentation is elementary and intuitive, and we omitted a number of 
technical details. 

The structure of our paper is as follows. The first part of this review is devoted 
to the description of the Yang-Lee theory for equilibrium phase transitions. A first 
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paragraph is devoted to the question of the role played by the thermodynamic 
limit. Then the general Yang-Lee theory for the grand-canonical partition function 

is discussed. The celebrated Circle Theorem is revisited, as well as its generalization 
to different models. Several aspects of the characterization of the critical behavior 
in terms of the Yang-Lee zeros are analyzed (density of zeros near criticality, Yang- 
Lee edge singularity, finite size effects). A brief discussion of the relation between 
the Pirogov-Sinai theory and the Yang-Lee zeros for first-order transitions is given. 
Then the problem of the Fisher zeros of the canonical partition function is revisited. 
Finally, the Potts model is discussed in order to compare several possible types of 
description of the phase transitions within the Yang-Lee formalism. 

The second part of this review is devoted to the question of extending Yang-Lee 
formalism to NESS phase transitions. In a first paragraph, the description in terms of 
a master equation is given, and the question of detailed balance discussed. Then the 
problem of a possible candidate for nonequilibrium steady-state partition function 
is approached. Applications to driven-diffusive systems, reaction-diffusion systems, 
directed percolation models, and systems exhibiting self-organized criticality are 
reviewed. The connection with equilibrium system with long-range interactions is 
also considered. Finally, some conclusions and perspectives are given in the last 
paragraph. 

PART I: EQUILIBRIUM PHASE TRANSITIONS 



2. Phase transitions and the thermodynamic limit 

When one varies a control parameter (for example, the temperature or an external 

field), a physical system can exhibit a qualitative change of its equilibrium state, 
i.e., it can undergo a phase transition. In the framework of thermodynamics such a 
phase transition shows up as a discontinuity or a singularity of physical observables 
(e.g., specific heat, susceptibility, etc.) as functions of the control parameter. A 
characteristic thermodynamic potential or some of its derivatives are discontinuous 
or singular, thus nonanalytic at the transition point. 

How is it possible to understand this phenomenon from a microscopic perspec- 
tive? From the point of view of the equilibrium statistical mechanics, the thermody- 
namic potentials can be expressed in terms of some properties of the microstates in 
the phase space of the system, that are determined by the microscopic interactions 
between the constituents of the system (particles, spins, etc.), as well as by the 
externally imposed constraints like a fixed value of a control parameter. More ex- 
plicitely, the characteristic potential is proportional to the logarithm of the partition 
function of the system. This partition function is a sum of the statistical weights 
over different configurations that are accessible to the system in the phase space 
under the given constraints, and these statistical weights are positively defined an- 
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alytic functions of the control parameter. If the sum contains a finite number of 
such statistical weights, then the partition function is finite, and its logarithm, i.e., 

the corresponding characteristic potential is an analytic function of the control pa- 
rameter, and therefore no phase transition is possible. To obtain a nonanalyticity of 
the characteristic potential, the partition function should thus necessarily become 
zero for a certain value of the control parameter. However, as justified above, this 
cannot happen in a finite system. A nonanalytic behaviour could only be obtained 
for a system containing an infinite number N of constituents, and thus one has to 
work in the so-called thermodynamic limit, in which both N ^ oo and the volume 
of the system V — > oo, and such that the density n = N/V is kept constant. 

One may object that all the systems in nature are finite. However, because 
of the huge number of its components (of the order of Avogadro's number), any 
macroscopic system behaves, from an experimental point of view, like an infinite 
system (except eventually in a very narrow domain in the vicinity of the transition 
point). The effects of the finite-size of the system on the appearance and properties 
of the phase transition will be subject of further discussion in Sec. 10. 

To summarize, the presence of a phase transition, from a statistical mechanics 
point of view, should be related to the vanishing of the partition function for a 
certain value of the control parameter. Thus one has to look for the zeros of the 
partition function, and such zeros should only show up in the thermodynamic limit. 

Let us briefly comment now on the properties of the interaction potential be- 
tween the constituents of the system that allow for the existence of a correct ther- 
modynamic limit. ^ Roughly speaking (see Refs. 1, 2, 3, 4, 5 for further technical 
details), for systems with two-body central interactions, the interaction potential 
u{r) (with r the distance between the particles) has to obey the following three 
conditions: 

(i) The intermolecular forces have to approach zero rapidly enough for large inter- 



''For example, consider a system in equilibrium with a heat bath at temperature T (that repre- 
sents a control parameter in this case). The free energy (which is the characteristic potential) is 
proportional to the logarithm of the canonical partition function, that is defined as a sum over 
all the microstates a of the Boltzmann weight exp(— Ba/fes^'), where Ea is the energy of the 
microstate a. 

''For such a system, under some supplementary constraints on the regularity of its frontier, taking 
the thermodynamic limit A'' oo, V ^ oo at fixed n = N/V implies that: (i) The various sta- 
tistical ensembles are equivalent, and the relevant thermodynamic parameters are thus uniquely 
defined. For example, the pressure has the same thermodynamic limit, positively definite, in both 
the canonical and grand-canonical ensembles, (ii) The values of the thermodynamic potentials (en- 
ergy, free energy, grand-canonical potential, enthalpy, etc.) per particle are finite. The system enjoys 
the additivity property, which means that its total Hamiltonian is the sum of the Hamiltonians 
of any combination of macroscopic parts taken separately (i.e. , one can safely neglect the surface 
interaction energy between macroscopic component parts) - and so do all the thermodynamic 
potentials (free energy, enthalpy, etc.). (iii) The thermodynamic stability of the system is ensured, 
i.e., the thermodynamic potentials enjoy the required conditions of convexity/concavity. For ex- 
ample, the compressibility coefficient is positive, which means the pressure is a non-decreasing 
function of the particle-number density n. (iv) The effects of the boundaries (e.g., surface energy 
or tension) are extinguished, and we speak only of bulk effects. 
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particle separations, i.e., in terms of the interaction potential, |w(r)| ^ Ci/r'^'^'^ as 
r — > oo (with Ci > and e positive constants), where d is the dimension of the 
system. Such systems arc currently called short-range interaction systems. 

(ii) u{r) has to have a repulsive part increasing sufficiently rapid at small inter- 
particle distances (preventing the system from collapse at high particle number 
densities). Throughout this paper, we shall consider systems with a hard-core in- 
terparticle repulsion at small distances, u{r) = oo for r < ro ~ 6^/'' (where b is 
the single-particle excluded volume) , which means that a system of volume V can 
accommodate a finite maximum number of particles M = V/b.'^ 

(iii) Finally, the interaction potential has to be everywhere bounded from below, 
u{r) ^ —uq whatever r (with Uq a positive constant). 

Note a rather general result concerning the one- dimensional equilibrium sys- 
tems with short range interactions, namely van Hove's theorem according to 
which no phase transition is possible in such systems, in contrast to what happens 
generically in corresponding nonequilibrium systems, see Ref. 8 for a brief review. 
A recent critical discussion ^ of this result allowed to highlight some exceptions to 
this theorem, but we shall not be concerned here with these rather pathological 
situations. 

Many physical situations can be modeled by interaction potentials with the 
characteristics (i) - (iii) above. There are, however, a few exceptions. The most 
salient example is that of the systems with long-range interactions (or non-additive 
systems), which include, e.g., systems with gravitational or Coulombian forces (see 
Ref. 6 for a modern review of their (thermo) dynamical properties studies). Such 
systems exhibit an unequivalence of the statistical ensembles, and it is not clear 
yet how to apply the concepts of the Yang-Lee theory (as described below) to their 
phase transitions. We will not address them further here. A limiting case of such 
systems is that of mean-field interactions, on which we will briefly comment in 
Sees. 7.2, 12, 13.3. 

3. Yang-Lee zeros of the grand-canonical pcirtition function: the 
general framework 

We shall address here two problems, namely: 

(A) The location of the phase transition point: How, when knowing the 
partition function (in the thermodynamic limit), can one locate a phase transition 
point by investigating the zeros of this partition function with respect to the control 
parameter of the system. 

(B) The characteristics of the transition: How can one extract information on 



'^As shown in Ref. 3, one can relax the hard-core condition to it(r) ^ C2/r'^~^^ as r — > (with C2 
and £ positive constants), but we will not consider this situation here. 
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the nature of the phase transition (e.g., if it is a discontinuous or a continuous one) 
from the properties and distribution of these zeros. 

The presentation in this section follows, with slight modifications, the same gen- 
eral lines of thought as in Ref. 10. 

(A) The location of the phase transition point. Let us start with this 
first problem and concentrate on the case of a generic system with short-range 
interactions as described above, for which the thermodynamic limit is well-defined. 
For concrctencss, and following the original papers of Yang and Lec ^^i^^ consider 
here the grand- canonical description of the system. The corresponding partition 
function, for a given volume V and a fixed temperature T, is expressed as: 

M 

Ev{T,z)=J2Zm{T)z^, (1) 

m=0 

where M = V/b (with b the single-particle excluded volume) is the maximum num- 
ber of particles that can be accommodated in the system (as determined by the 
hard-core part of the binary interaction potential); Zm{T) is the canonical parti- 
tion function of the system with fixed number of particles m; and finally z is the 
fugacity of the system, 

z = exp(/i/fcBT) , (2) 

which is obviously a real, positive quantity expressed in terms of fi, the chemical 
potential of the system in contact with an external reservoir of particles; fcs is 
Boltzmann's constant. 

One notices that Sy(T, z) is a polynomial of M-th order in the fugacity z. Let us 
consider its roots, i.e., the solutions of Ev{T,z) = 0; according to the fundamental 
theorem of algebra, there are M such roots, Zi = Zi{T), with i = 1, M. Moreover, 
because the coefficients of all the powers of z in the expression of Sy (T, z) are real 
and positive, all these roots Zi{T) appear in complex- conjugated pairs in the complex- 
fugacity plane, away from the real, positive semi-axis. One can express Hy(T, 2) in 
terms of these roots, 



Ev{T,z)=Cl[ 



M 

7 

(3) 



i=l 



ZiiT) 



(where ^ is a multiplicative constant, that we shall ignore in the foregoing), and the 
corresponding finite-size grand-canonical potential fiy = — /cB^lnSy = —Py V 
leads to the following expression of the finite-size pressure Py: 



1 ^ 

Pv{z) = kBT-J2^n 



M 

7. 

(4) 



^i(r). 



I 



I 
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Note that throughout this Section we shaU consider the temperature T as a fixed 
parameter, and therefore aU the quantities hke Pv{z), p{z), P{z), ^{z), il>{z), A(s), 
as well as C (sec below) are temperature- dependent. 

Let us consider now the complex extension of Pv{z), that is defined through 
Eq. (4) above for all complex z with the exception of the points Zi. Prom the standard 
theory of the complex- variable functions, one can deduce that Pv{z) is analytical 
(infinitely differentiable) over any region of the complex- 2; plane that is free of zeros 
of the partition function. Therefore, a nonanaliticity of Pv{z) at a complex point 
zq appears if and only if in any arbitrarily small region around zq one finds at 
least one root of the partition function (or, in other terms, if and only if there 
is an accumulation of the roots of the partition function in the vicinity of Zq). If, 
moreover, zq lies on the physically accessible real positive semi-axis, this corresponds 
to a phase transition in the system. 

Once more, one realizes that such conditions cannot bo accomplished in a finite- 
size system. Let us then turn to the thermodynamic limit, where one might even- 
tually expect a possible accumulation of the roots of the partition function towards 
the real positive semi-axis. 

As V and therefore M = V/b increase, the location of the roots Zi changes, and 
in the thermodynamic limit they accumulate in a certain region C of the complex- 
z plane, with a local density p{z). Of course, in view of its significance, p{z) is a 
real-valued, non-negatively defined function on C, identically zero for any z outside 
C, which is normalized as 



and thus it is integrable over any bounded region of the complex-z plane. ^ Moreover, 
the symmetry property of the roots with respect to the real axis is preserved in the 
thermodynamic limit, and leads to the following property of the density of zeros: 



(where (...)* denotes complex-conjugation), i.e., C is symmetric with respect to the 
real axis. Note also that the region C depends on the temperature T, and, of course, 
its characteristics are specific to each system. For all the points z outside the region 
C, one can define the thermodynamic limit of the complex extension of the pressure, 
P{z) = limy_>ooPy(-2), a complex-valued quantity expressed as 





p{z) = p{z*) 



(6) 





"^As already mentioned in Sec. 2, we shall consider here only systems with a hard-eore interparticle 
repulsion at small distances, i.e., with a non-zero single-particle excluded volume 6. Most of the 
litterature on the subject is limited to this situation; see Ref. 22 for an exception. 
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which is, of course, a multi-valued function. Let us consider its real part (up to a 
factor ksT) 



^{z) = ^HeP{z) = J dz'p{z')\n 



'-5 



(8) 



Since In |z| (with z = x + iy; x, y E R) is the Green function of the two-dimensional 
Laplacian A = d"^ /dx^ + d"^ jdy^, by applying the Laplacian to the above equation, 
one finds that: 

p{z) = i- A^(^) . (9) 

This means that <^{z) is the electrostatic potential associated to a distribution of 
charges of density p{z). This analogy is very useful, since it allows the direct tran- 
scription of well-known results from electrostatics. In particular, since p{z) is in- 
tegrable even on regions containing parts of C, the function ^(z) can be extended 
through continuity over the entire complex-^ plane. The equipotential surfaces are 
thus given by ^5(2)= const ant. Still in this analogy with electrostatics, the imaginary 
part of P{z)/kBT, 

V'(^) = j^^^P^^) (10) 

(defined modido 27r) determines, through the condition '(/'(z)=constant, the lines of 
force of the electrostatic field generated by the distribution of charges p{z), and the 
intensity of the field is given as \J^{z) and is, of course, discontinuous on C, see 
point (ii) below for further details. Note also that 

^p{z) = ip{z*) and ij{z) = -tp{z*) , (11) 

which result directly from the symmetry property (6) of the distribution of zeros. 

Suppose now that the known partition function of the system indicates us that 
ip{z) has, for example, two distinct analytical expressions ipi{z) and (p2{z) in two 
different regions of the complex-z plane. But since (p{z) has to be continuous over 
the entire plane, there should be a matching region between these expressions, and, 
of course, this region can be nothing else than C. Indeed, since (pi^2 arc different 
functions, there appears the possibility of a nonanaliticity of (p (and thus of P{z)) 
at the matching points. Then the location of C is given by the condition: 

ipi{z) = ip2{z) ,i.e., nePi{z) = 7^eP2(^) • (12) 
c c c c 

If the region C intersects the real positive semi- axis at some point zo, and/or 
there is an accumulation of zeros in the vicinity of Zq, then we have a phase transition 
in the system, and 

7^ePl(zo) =7^eP2(zo), (13) 



the real part of the complex-valued pressure is continuous at the transition point. 
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Of course, depending on the shape of C, several transition points may appear, 
for different values of the fugacity; one can also encounter the phenomenon of mul- 
tiple phase coexistence - like, for example, the appearance of triple points through 
the coalscence, at some temperature, of two simple (two-phase) transition points. 
Moreover, recall that C depends on the temperature, and from here there appears 
the possibility of the existence of a critical temperature Tc, such that for any T < Tc 
the zeros of the partition function accumulate in the vicinity of the real positive 
semi-axis (i.e., there is a phase transition in the system), while for any T > Tc the 
zeros of the partition function are no longer accumulating in the vicinity of the real 
positive semi-axis, i.e., a phase transition is no longer possible. 

We thus answered the point (A) concerning the formal frame for the description 
of the appearance of the phase transitions - at least, for the time being, in the 
grand-canonical ensemble. 

(B) The characteristics of the transition. In order to answer the second 
question, namely the quantitative connection between the nature of the phase tran- 
sition and the distribution of zeros of the partition function, let us suppose for the 
moment that C is a smooth curve in the vicinity of the transition point zq. As it will 
be seen below, on some concrete prototypical examples of physical systems, this is 
rather often the case for the zeros of the grand-canonical partition function in the 
complex-fugacity plane - the so-called Yang-Lee zeros. 

Following the arguments in Refs. 10, 13, 14, 15, 17, 18, 19, let us consider a 
parametrization of the curve C in the vicinity of zq, with the parameter s measuring 
the anti-clockwise oriented distance from zq along the curve (s = at z = zq), see 
Fig. 1. 



A y = Jmz 



C 



t 



ipiiz) -f ii'i{z) 




n 



X — Tie z 



Fig. 1. Schematic 
transition point zq, 



representation of the location C of the Yang-Lee zeros in the vicinity of a 
meant to illustrate the notations in the main text. 
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Let A(s) denote the line density of zeros along C. Then, according to Gauss' 
theorem, one has for a point z on the curve 



^ [V<p2(2;) - V<^i(z)] • n 



= A(s) , (14) 

c 



where h = n{s) is the unit vector normal to C (oriented from region "1" to region 
"2") at that point z. In view of the Cauchy-Riemann properties of the analytical 
functions Pi,2(^), one has 



where 



V<pi2(2)-n = ViJio{z)-i (15) 
c c 



^1,2 (^) = T:Kf;IrnPi,2{z) (16) 



are the imaginary parts of P\^2{z)/kBT (defined modulo 2it), and i = i{s) is the 
unit vector tangent to C. This leads finally to 



c 
(17) 



where d/ds denotes the directional derivative along the curve. In particular, at the 
transition point zo, the discontinuity in the directional derivative of the imaginary 
part of the complex-valued pressure is determined by the local linear density of 
zeros. 



^ {ImP2 (z) - ImPi (z) 



27rfcsT ds 



m ■ (18) 



The relationships (12), (13), (17), and (18) will allow us to establish the nature 
of the phase transition, i.e., whether it is discontinuous (or first-order) or continuous 
- i.e., second- or higher-order in the classical Ehrenfest classification scheme. 

Let us consider the Taylor expansion around zq of the complex-valued pressure 
on both sides of the curve C, i.e., 

^ -PiA^) = j^Pizo) + ai,2{z - Zo) + h,2{z - ^o)' + Oiiz - zo)^) , (19) 



In order for the pressure to be real on the real z axis, all the coefficients of the 
development have to be real. From the condition (12), one finds that the equation 
for the curve C is given by 

(as -ai){x- zo) + (62 - 61) [(a; - ^0)^ - y'] + He [0{z - zo)^] = , (20) 

where x and y are the real, respectively the imaginary part of z. Several situations 

may appear. 

(i) First-order phase transition. If 02 ^ ai then the complex-valued pressure has 
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a discontinuity in its first derivative and the transition is first-order. If, moreover, 
&2 7^ &i, then in the vicinity of the transition point zq the curve C is a hyperbola, 

y^ = (x-Zof + ^P^ix-Zo) (21) 

02 - h 

whose tangent in zq is parallel to the imaginary axis, i.e., C crosses the real axis 
smoothly at an angle ■it/2. Using Eq. (18) we find the density of zeros at zq, 

A(0) = ^, (22) 

i.e., the density of zeros at the transition point of a first-order phase transition is 

nonzero. 

(ii) Second-order phase transition. If 02 = ai, but 62 5^ 61, then the curve C obeys 
the equation 

y = ±ix-zo), (23) 

i.e., in the vicinity of zq it consists of two straight lines that make an angle of ±7r/4 
with the real axis (and 7r/2 between them) and meet at zq. Prom Eq. (17) we find 
that 

X{s)=^-^^\s\+0{s^), (24) 

TT 

i.e., the desity of zeros is decreasing linearly to zero when approaching the transition 
point zo {s — 0). 

(iii) Higher-order phase transitions. If the discontinuities appear at higher order in 
the derivatives of the complex-valued pressure, then one can repeat the above type 
of reasoning to find the equation of C and the density of zeros in the vicinity of the 
transition point zq. In general, if the transition is of n-th order (n ^ 3), then the 
density of zeros is zero at the transition point, A(s) ~ and the curve C does 
not cross smoothly the real axis, but approaches it at an angle ±7r/2n from above 
and below. 

Resuming, we managed therefore to respond the two main questions: 

(A) The accumulation of zeros of the partition function along the (physically aces- 
siblc) real, positive semi-axis of the complexified fugacity z indicates the location 
of the phase transition point(s); 

(B) The density of zeros near such an accumulation point determines the order of 
the transition (according to Ehrenfest's classification scheme) at that point. 

In this section wc discussed the zeros Yang-Lcc zeros of the grand-canonical 
partition function in the complex-fugacity plane. One can, of course, address the 
following legitimate question: how can one extend this type of formalism in order 
to study phase transitions in other statistical ensembles (e.g., the canonical one). 
Indeed, in view of the equivalence of these ensembles in tlic^ thermodynamic limit, 
one is entitled to expect similar results on the location and characteristics of the 
phase transitions whatever the ensemble used in the description of the system. 
Before addressing this important problem in Sec. 12 below, let us add a few more 
relevant comments on the Yang-Lee zeros. 
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4. Grand-canonical partition function zeros and the equation of 
state 

We saw above that knowing the distribution p{z) of the zeros of the grand-canonical 
partition function in the thermodynamic limit, one can compute the complex-valued 
pressure P{z) through Eq. (7), 

P{z) = kBTxiz) , (25) 

with 

X{z)^ I dz'p{z')ln{l-^) , (26) 

that is analytic everywhere outside C. Also, the complex- valued particle- number 
density n{z) follows from the complex extension of 



niz) = z'4^ = z fdz'^, (28) 



as: 



dz 

for all complex z outside C. When z is real and positive (and outside the possible in- 
tersections of C with the real positive semi-axis that corresponds to phase-transition 
points), Eqs. (25) and (28) (with x(z) given by Eq. (26) as a functional of p{z)) 
are parametric expressions of the equation of state of the system. Recall also that 
p{z) and C, and therefore x{z) are temperature-dependent, see Sec. 3, and thus the 
structure of these parametric equations is actually quite intricate. Of course, one 
can, in principle, eliminate z between the two equations and find the explicit form 
of P = P{n, T) (for the different regions of the real-positive semi-axis of z outside 
C, i.e., free of transition points). 

However, one has to realize that the problem of determining the distribution of 
complex zeros p{z) is an extremely diSicult task, even in the simplest known-cases 
like, e.g., a one-dimensional lattice gas with nearest-neighbour attractive interac- 
tions ^'^ , a gas of hard rods a mean-field lattice gas ^o,4i_ jg q£ interest, therefore, 
to address the following problem. Suppose that we are given an certain equation of 
state P = P{n,T); using the parametrization of Eqs. (25) and (28), one can find 
a closed nonlinear first-order difii'erential equation for the function x{z), which, by 
integration, leads generically to a functional equation of the form 

Hx,z)=0, (29) 

(for real positive z). Let us extend it, by definition, to the whole complex-z plane. 
Would it then be possible to solve it and, once x(^) is found, to determine the region 
C of accumulation of zeros of the grand-canonical partition function? Supposing that 
C is a smooth curve, would it be then possible to find, through Eq. (26) or Eq. (17), 
the corresponding density of zeros p{z) of the grand-canonical partition function of 
the system? 
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This is what is currently called the inverse problem (in analogy with the inverse 
problem in electrostatics), and it has been addressed in several papers on various 
specific systems, sec Rcfs. 20, 21 (two models of Ising fcrromagncts) , 22 (a Tonks 
gas of hard rods and a gas with a weak long-range repulsion), 31 (Tonks gas), 23, 24, 
27, 28, 29, 30 (the van der Waals gas), 25 (a lattice gas with a hard-core repulsion 
extended over several lattice sites), 32, 33 (Takahashi lattice gas), 34 (point particles 
with logarithmic; interactions), and even in an experimental setup, see Ref. 35 (for 
a two-dimensional Ising ferromagnet). 

The answer to this issue is related to the considerations in the previous section 
and, being rather technical, will not be given here in detail; we refer the reader 
to the set of four papers by Ikeda, Refs. 36, 37, 38, 39, for a rigorous approach 
of the problem (and also some illuminating examples, like various systems with 
circular distributions of zeros, see Sec. 7 below; the ideal Fermi-Dirac gas; and the 
ideal Bose-Einstein gas). Roughly, an equation of type (29) has as solutions one or 
several complex- valued functions, which may have one or several Riemann surfaces 
(i.e., are in general multi- valued functions), with branching points corresponding to 
dz/dx — 0. One has to try to match these different Riemann sheets, along the curve 
C, so that several mathematically and physically necessary conditions are fulfilled 
by the resulting patch- function xi^)- (i) The branching points belong to C; (ii) The 
real part of x(z) {ip{z) in the notations of the previous section) has to be continuos 
in the whole z plane; (iii) xi^) is real and continuous on the real positive semi-axis 
and has the value given by Eq. (25); (iv) The integration constant result from the 
infinite-dilution (ideal gas) limit, x{^) ^ -2 for ^ — > 0; (v) The quantity p{z), as 
defined by Eq. (17), has to be real- valued and non- negative for all z on C. The 
discontinuities in the imaginary part of x{z) {i^iz) in the notations of the previous 
section) appear thus across the charged contour C originating from matching the 
different Riemann sheets of the multi- valued solutions of (29). In the particular case 
when the solution of (17) is a single multi- valued function, the domain C of the zeros 
of the partition function reduces to the branching points. One has to note, however, 
that for a general functional relation (29) there is no guarantee of the uniqueness of 
the domain C, and also that rather often the explicit construction of C, and thus the 
calculation of p{z) are not possible, except for limiting cases, e.g., in the vicinity of 
the branching points. 

5. The lattice gas and the Ising model in a magnetic field 

Before proceeding in the next sections with a discussion of the domain C of accumu- 
lation of zeros of the grand-canonical partition function, let us recall a well-known 
classical result of equilibrium statistical mechanics, namely that the problem of a 
discrete lattice gas (with nearest-neighbour interactions) is mathematically equiv- 
alent to the problem of an Ising model (with nearest-neighbour interactions) in an 
uniform magnetic field, see, e.g., Refs. 12, 4 for a detailed description. In partic- 
ular, the expression of the grand- canonical partition function of the lattice gas is 
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identical to the expression of the canonical partition function of the king model in 
magnetic field (at the same fixed temperature T). The search of the complex roots 
of the grand-canonical partition function of the gas in the complcx-fugacity plane 
amounts to the search of the complex roots of the canonical partition function of 
the Ising model in the plane of the complex variable z = exp{—2H/kBT) (in the 
appropriate units for the field and the temperature), that is related to a complex 
magnetic field intensity H. 

Therefore, based on this one-to-one correspondence, in the foregoing we will 
speak either of a discrete lattice gas, or of an Ising model in a magnetic field. 

6. Yang-Lee zeros and the transfer matrix formalism 

For one-dimensional Ising lattice systems with finite-range interactions, in a mag- 
netic field it is customary to compute the partition function using the transfer- 
m,atrix technique. The canonical partition function for A'^ spins can be written as: 



where Aj, with i = 1, fc, are the k eigenvalues of the k x k transfer matrix of the 
system. Besides being dependent on the reduced coupling constants (i.e., the cou- 
pling constants divided by fcsT) between the spins, these eigenvalues depend on the 
value of the magnetic field H , or, cquivalcntly, on the fugacity z = exp(— 2i7/fcBT), 
Aj = \{z), i = 1, k. Let us consider the free energy per spin, F = —kBTlnZ]\j /N , 
and its extension to the plane of the complex fugacity. In the thermodynamic limit, 
only the eigenvalue of highest modulus contribute to the free energy. Suppose now 
that in two different regions of the complex- ^; plane two different eigenvalues, Ai(2:), 
respectively A2 (2) assume the largest modulus. Since the real part of the free energy 
(per spin) has to be continuous throughout the whole z plane - see the discussion 
in Sec. 3 in the light of the correspondence between Ising systems and lattice gases. 
Sec. 5 - it follows that the location C of the zeros of the partition function in the 
complex- 2; plane is given by the condition of matching of the modulus of these two 
eigenvalues of the transfer matrix: 



which is the equivalent of Eq. (12) in the present frame. Of course, in view of 
the van Hove theorem (Sec. 2) for the systems with short-range interactions that 
we arc considering here, no phase transition is possible, i.e., C does not have any 
accumulation point on the real positive semi-axis for such systems. This matching 
condition between eigenvalues might be useful for the construction of C, see, e.g., 
Ref. 45, 110 for a few examples. A more rigorous discussion of this type of approach 
for a wide class of lattice gases can be found in Ref. 48. 



Zj, = X^ +... + X- 



N 
'k ' 



(30) 



|Ai(2)| = \X2{z)\ 
c c 



(31) 
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7. The Circle Theorem for the Yang-Lee zeros 

The general frame discussed in Section 3 above is applicable to all the systems - 
within, of course, the mentioned restrictions that refer, essentially, to the existence 
of the thermodynamic limit, sec Sec. 2. However, as already mentioned, the shape of 
the domain C of the accumulation of zeros of the grand-canonical partition function 
in the complex-fugacity plane, besides being temperature-dependent, is specific to 
each system, and a priori one might expect it to have a rather intricate structure. 
It is therefore a rather amazing result that, for a quite general class of systems, this 
region C can be proven rigorously to lie on the unit circle \z\ = 1 in the complex- 
fugacity plane. This is the celebrated Circle Theorem, that was first demonstated 
by Yang and Lee, see Rcf 12. for a rather restricted class of Ising systems. The 
beauty and simplicity of this result attracted a lot of research afterwards, and led 
to its generalization to many other situations and models, see below. Note, however, 
that no general statement can be made about the density of zeros on the unit circle, 
which is characteristic to each system and, of course, temperature-dependent. 

7.1. The original Yang- Lee Circle Theorem 

The system originally considered by Yang and Lee (see Refs. 12, 72, 26) consists of 
an ensemble of N Ising 1/2-spins at {i = 1, N): 

(i) that arc placed on a d-dimcnsional lattice, 

(ii) with pair ferromagnetic interactions between them, of coupling constants Jij ^ 



leads to a canonical partition function = Zj\[{zi,Z2, zn) that is a multinomial 
in the fugacities 



{Hi being the value of the magnetic field acting on the spin cr,). In the particular 

case of a homogeneous magnetic field Hi = H for all i, the partition function is 
simply a polynomial of degree A'' in the fugacity z = exp{—2H/kBT). Under the 
supplementary hypothesis 
(iv) \zi\ ^ 1 for all i = 1,...,N 

one can prove the Circle Theorem: the zeros of the partition function in the 2N- 
dimensional space of complex Zi lie all on the unit circle. 





(32) 



i<j 



Zi = exp{-2Hi/kBT) 



(33) 



Zn = 0^ \zi\ = \Z2\ = ... = |-2jv| = 1 , 



(34) 



°Of course, alternatively, one can consider the equivalent lattice gas model, which corresponds to 
an uniform magnetic field. 
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which amounts at saying that these roots correspond to strictly imaginary (or zero) 
values of the magnetic field intensity Hi, i = 1,...,N. In the particular case of a 
homogeneous magnetic field, all the A'' roots of = lie on the unit circle = 1 
in the complex-z plane, i.e., appear for N strictly imaginary (or zero) values of 
H. The theorem remains valid in the thermodynamic limit, when the zeros form 
a dense set C on the unit circle, characterized by the linear density A. Of course, 
an immediate consequence of this theorem is the well-known classical result that 
a phase transition might appear in such a system, at a fixed finite temperature T, 
only in zero magnetic field. 
Note that: 

(i) The dimensionality of the lattice does not play any role in the demonstration 
of the theorem, which is thus valid in any dimension. Moreover, the size and the 
structure of the lattice, its regularity, periodicity, translational symmetry do not 
play any role either in the obtention of the result. 

(ii) The ferromagnetic interactions are not restricted to first neighbours. However, 
they should decrease sufficiently rapidly with the distance in order to ensure the 
existence of the appropriate thermodynamic limit, see the discussion in Sec. 2. 

Indeed, the demonstration of this theorem that will not be reproduced here, 
rely on some properties of the Hamiltonian that are independent of the dimension 
and of the range of the interactions; besides the ferromagnetic character of the 
interspin interaction, an essential necessary ingredient is the spin-reversal symmetry 
of the system, i.e., the invariance of its Hamiltonian with respect to global inversion 
of the spins and of the (inhomogeneous) magnetic field, which leads to the following 
symmetry property of the canonic partition function: 

1- (35) 

Later on Asano, in Refs. 53, 54, introduced a new, more general technique for the 
study of the location of the zeros of the grand-canonical partition function. This 
method was often used to prove the Circle Theorem for various systems (see below), 
but it was also extended by Ruelle, see Ref. 73, 74, so as to permit statements about 
regions other than the unit circle. 

7.2. Generalizations of the Circle Theorem to other systems and 
models 

We shall briefly present below some of the systems for which it was shown, using 
various methods, that the Circle Theorem is valid. 

(A) Modified Ising ferromagnets, which include: 

(i) Ising ferromagnets of arbitrary spin, see Refs. 49, 52, 50. The basic idea for the 
extension is that the higher-spin "atom" is a cluster of 1 /2-spin "atoms" , that are 
coupled together through a suitable ferromagnetic interaction. This representation 
leads to an effective, temperature-dependent, 1/2-spin Ising Hamiltonian and from 
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here to the Circle Theorem. Once demonstrated this result for any finite value of 
the spins, one can extend it to the classical limit or infinite-spin limit, see Ref. 59. 

(ii) Diluted Ising ferromagnets. Consider a lattice whose sites are either occupied 
by a spin 1/2 (with a probability p), or by a non- magnetic atom (with probability 
1 — p); p is thus equal to the concentration of spins on the (infinite) lattice. As shown 
in Refs. 60, 61, for a ferromagnetic coupling between the spins, and in a homogenous 
magnetic field H, the zeros of the partition function correspond to purely imaginary 
(or zero) values of H. This result can be extended to higher-spin cases. 

(iii) Ising models with many-spin interactions. Ising lattice systems with interac- 
tions involving three or more sites - i.e., whose Hamiltonian comprises interaction 
terms of the form —Jij...kO'iCrj...ak, where i ^ j ^ ... ^ k correspond to different 
lattice sites, and Jij...k is the coupling constant between these sites - have been 
used to model a variety of physical situations, like, for example, some binary al- 
loys , lipid bilayers , and gauge-field theory models These systems have rich 
phase diagrams, and Monte-Carlo results indicate the presence of phase transitions 
at nonzero values of the magnetic field for a number of these systems, see, e.g., 
Refs. 87, 88. However, in view of their complexity, very few analytically rigorous 
results are known. In Ref. 55 it is shown that for a set of many-spin interactions of 
finite range, under some conditions on the coupling constants, and also some spin- 
inversion symmetry conditions, the Circle Theorem holds (at least) up to a certain 
temperature (that is determined by the coupling constant and does not depend on 
the number of spins in the system). The example of four-spin interactions shows 
explicitly that the zeros leave the unit circle at large-enough temperatures. Qual- 
itative arguments suggest that this is the case for other systems with multi-spin 
interactions, see also Ref. 83. 

(iv) Ising models on hierarchical lattices. In Ref. 89 a rescaling formalism is used 
in order to study the location of the Yang-Lee zeros for 1 /2-spin and 1-spin Ising 
chains on two- and three-dimensional Sierpinski gaskets, as well as on 3-simplex and 
4-simplex lattices. The Circle Theorem is shown to hold for nearest-neighbour fer- 
romagnetic interactions, but is no longer valid at high temperatures when four-spin 
interactions are included, see also point (iii) above. References 90, 91 address the 
problem of the Yang-Lee zeros of the Ising nearest-neighbour ferromagnetic model 
on a Cayley tree. The Julia set of the renormalization transformation of the model 
gives the thermodynamic limit of the partition function zeros distribution. This one 
lies on the unit circle, and is multi-fractal in nature, with a tempcratiire-dependent 
generalized dimension; below Tg the Julia set is the unit circle (i.e., has dimension 
D = 1), and above Tc the Julia set is a Cantor set (with a dimension D{T) < 1 
decreasing with the temperature) on the unit circle. When multi-site interactions 
are included, these zeros leave the imit circle, see Ref. 92. 

(v) Aperiodic Ising models. References 93, 94, 95 discuss various one- and two- 
dimensional examples of Ising spins with nearest-neighbour interactions on aperi- 
odic lattices. The interaction constant can take only two values, Ja.t, distributed on 
the lattice according to some generating rule. In ID case, for example, the actual 



2, 2008 17:43 WSPC/INSTRUCTION FILE Yang-Lee 



18 /. Bena, M. Droz, A. Lipowski 

distributions of the two coupling constants Jo, 6 along the chain is determined by 
an infinite "word" in the letters "a" and "b" which is obtained as the unique limit 
of certain two- letter substitution rules (like, e.g., the Fibonacci or the Thuc-Morse 
rules). If both Ja^b are positive (ferromagnetic coupling), then, although the zeros 
still lie on the unit circle, their distribution has a fractal structure (that, of course, 
disappears when Ja = Jb), as shown by the integrated density of the zeros along 
the circle. Thus such systems may be regarded as intermediate between regular and 
hierachical models. 

(vi) Mean-field Ising model. In Refs. 40, 41, 42, 43, is addressed the problem of the 

Yang-Lee zeros for the Temperley-Husimi model of a lattice gas, or for the equiva- 
lent mean-field Ising model in a magnetic field, which, of course, is known to exhibit 
a first-order phase transition, with a critical temperature T^. The partition function 
zeros lie on the unit circle in the complex-fugacity plane, and one can compute 
analytically the corresponding density. A direct comparison with Mayer's cluster 
expansion theory 9^,98,99 gj-^ows that the non-analiticity point in Mayer's series does 
not coincide with the true transition point as obtained from the Yang-Lee theory, 
see also related comments in, e.g., Refs. 11, 20, 100. See also Rcf. 101 for an ap- 
proach using the Yang-Lee zeros of a restricted partition function for the study of 
the metastable states of long-range (and their mean-field limit) Ising models. 

(B) Heisenberg spin models. The technique of studying the Yang-Lee zeros 
that was introduced by Asano in Refs. 53, 54 was also used by Suzuki and Fisher 
in Ref. 55, their studies being devoted to the extension of the Circle Theorem to 

quantum 1/2-spin systems in a magnetic field. More precisely, they considered the 
following fully anisotropic Heisenberg Hamilonian for a system of N spins of com- 
ponents crf'^'^, i = 1, ■■■,N, in a magnetic field: 

^ = - E (-^^^ + + 

-Y,{Hfaf+Hya« + H^a!), (36) 

i 

where J^j^'^ are the (anisotropic) coupling constants and H^'"^'^ are the (non- 
uniform) components of the magnetic field in the x,y and z directions. Consider 
that H^ = Hfii, with /ij non- negative constants (for all i = 1,...,N), and that the 
coupling constants obey the following set of inequalities 

Jij > I J^l and J^j > for all i ^ j = 1, N , (37) 

which correspond to the so-called ferromagnetic Heisenberg model with domi- 
nant z ~ z coupling. Then one can prove that the zeros of the partition function 
Zm — Tr [exp(— TY/i^sT)] as a function of H lie all on the imaginary- axis for all 
N, and do so also in the thermodynamic limit. This result was generalized further in 
Ref. 55 for quantum models with spin greater than 1/2 (note that some preliminary 
results have already been obtained in Ref. 56), and also to the classical Heisenberg 
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system of infinite spins, see Ref. 59; see also Ref. 61 for a temptative relaxation of 
the corresponding conditions (36) on the coupling constants in the classical limit. 
In Ref. 59 it was also shown that the Circle Theorem holds for a two-sublattice an- 
tiferromagnet of arbitrary spin, in terms of the staggered magnetic field intensity. 

(C) Monomer-dimer models. Consider a lattice of N sites or vertices, connected 
through N{N — l)/2 edges; one associates to the edge connecting the vertices i and 
j a weight Wij , and to a vertex i a weight Xi (of course, Wij and Xi are non-negative 
numbers for all i ^ j = 1, ...,N). A dimer- monomer covering of this lattice corre- 
sponds to (i) placing dimers on the edges such that no two dimcrs share a common 
vertex (i.e., no vertex has more than one dimer); (ii) placing monomers on all the 
remaining vertices, i.e., the vertices that are not adjacent to a dimer. The statistical 
weight of such a covering is then the product of all the weights Wij of the edges 
covered by the dimcrs times the weights Xi of all the vertices covered by a monomer. 
The partition function Zn = Zj\[{{wij} , {xi}) of the system is then the sum of the 
weights of all the possible coverings. These models correspond to a huge variety of 
physical situations, see, e.g., Ref. 64 for an overview. 

References 62, 63, 64, 61, 65 address the problem of the zeros of the partition 
function of the monomer-dimer system. It is shown that, for given non- negative 
Wij, Zn cannot be zero if TZe{xi) > for all i = 1, or if TZe{xi) < for all 
i = 1, A^. In particular, ii Xi = x for all i = 1, N, this leads to the conclusion 
that the zeros of Z^ appear on the imaginary- a; axis (or at x = 0), a result that 
holds also in the appropriate thermodynamic limit. This corresponds to the absence 
of a phase transition in the monomer-dimer system, except, eventually, at .x = 0, 
which corresponds to a maximum density of dimers. Connections with Ising and 
Heisenberg spin models are also discussed. 

(D) Ferroelectric models. In Refs. 55, 57, 58 it was shown that the Slater- 
type models for ferroelectricity can be mapped onto Ising models with at most 
four-spin interactions. Therefore, the Yang-Lee zeros lie on the unit circle at low 
temperatures, but leave it at the critical temperature Tc- A modified Slater model, 
and an antiferromagnetic model are also discussed from the perspective of the Yang- 
Lee theory, and their distributions of zeros on the unit circle, for finite-size lattices, 
are investigated numerically. 

(E) Quantum fields. References 66, 67, 68, 69, 70, 71 show that some Euclid- 
ian quantum fields can be approximated by generalized Ising models on a lattice, 
with a given probability distribution of the spins (that may be discrete or contin- 
ues, but has to obey certain symmetry properties), and with a self-interaction of 
the spins. In particular, the Yang-Lee Circle Theorem is shown to hold for the con- 
sidered systems, see the above-cited references for further details. 
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For all these systems for which the Circle Theorem holds, using the notations 
in Sec. 3 one can say that the curve C lies on the unit circle \z\ = 1, and can be 
parametrized by the trigonometric angle 9. Correspondingly, the density of zeros 
A = X{9) = X{—6) (in virtue of the symmetry property (6)), and it is nonzero on 
two arches of the unit circle defined by ^ |^| ^ ti", where > is a function of 
temperature (see also Sec. 9.2). The normalization condition (5) for the density of 
zeros reads simply 

2 / \{e)de = b. (38) 
ho 

The complex pressure (7) becomes 



P = kBT [ X{0)ln{z^ -2zcose + l)de, (39) 
Jeo 

and the complex particle-number density (28) is 

n = 2z r \{0) , V . d0. (40) 
Jg^ ^ ^ z2-2^cos6' + l ^ ' 

Note, however, that the density of zeros X{9) is known analytically in very few 
cases (e.g., the one-dimensional Ising model with nearest- neighbour interactions, see 
Ref. 12; the mean-field Ising model, Refs. 40, 41), and was computed numerically 
in few other cases (see Ref. 21 for such a calculation for the Ising ferromagnets 
on a two-dimensional square lattice, and on a three-dimensional diamond lattice, 
respectively), and thus the practical use of these relationships is rather restricted. 



7.3. Griffiths inequalities 

Note also a very important by-product of the Circle Theorem, namely the 

demonstration of various Griffiths-type inequalities, i.e., inequalities between spin- 
correlation functions, see Refs. 50, 53, 54, 51, 66, 67, 69, 70, 75, 76. The simplest 
to demonstrate are the Griffiths-Kelly-Sherman 77,78,79,80 (qKS) inequalities. They 
were first obtained using the Circle Theorem in Refs. 50, 51, for an Ising 1/2-spin 
ferromagnet in a magnetic field iJ ^ 0, as: 

(ctactb) > (the first-type GKS inequafity) (41) 

{o'ao'bo'c'^d) ^ {o'ao'b) {o'co'd) (the second-type GKS inequality) , (42) 

where <ja,b.c.d arc the products of the spins inside the (multi-site) regions 
A,B,C,D of the lattice, respectively; (...) denotes the mean over the statistical 
ensemble. These inequalities were generalized, on the basis of Yang-Lee Circle The- 
orem, to higher-spin systems^'', Heisenberg ferromagnets^^'^^, and systems with 
arbitrary even spin distributions®''. 

The Griffiths-Hurt-Sherman*i (GHS) inequality 



{oaoboc) ^ {oa){oboc) + {oc){oaob) + {(yB){ocOA) - 2{a a) {ctb) {(Tc) (43) 
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was demonstrated for various systems, using the Circle Theorem, in Refs. 66, 69, 70. 
It has several interesting implications, the concavity of the average magnetization 
as a function of H; the monotonicity of the correlation length as a function of the 
magnetic field; the well-known result on the absence of phase-transition for H > 0; 
some inequalities for the critical exponents; see, e.g., Refs. 82, 76. More complicated 
inequalities on the correlation functions, involving the Circle Theorem, are treated 
in Refs. 69, 75. 

8. Further results on the Yang-Lee zeros location 

As already mentioned, the papers by Asano, Refs. 53, 54, and later on their elegant 
generalization by Ruelle, Refs. 73, 74, lead to a criterion for the location of the 
Yang-Lee zeros in the complex-fugacity plane for a fairly large number of systems. 
This criterion, as well as ad-hoc methods adapted for the system under study (e.g., 
the transfer-matrix formalism for the one-dimensional lattice systems, see Sec. 6), 
and sometimes numerical results allowed to approach several systems and to draw 
some general conclusions. A first remark would be that usually the distribution of 
Yang-Lee zeros in the complex-fugacity plane is well-behaved, i.e., in the thermo- 
dynamic limit the zeros are distributed densely on an ensemble of smooth curves. 

(A) Ising systems with multi-spin interactions. A rather general result refers 
to the fact that the Circle Theorem does not hold at high temperatures, see points 
(iii) and (iv) for modified Ising models in the previous section. A rigorous study 
of the regions in the complex-magnetic plane that arc free of zeros of the partition 
function is carried in Ref. 83 for two very general classes of multisitc interaction sys- 
tems, namely (i) systems having ferromagnetic interactions involving even numbers 
of sites, and (ii) systems with, again, ferromagnetic interactions involving even num- 
bers of sites and, in addition, either ferromagnetic or antiferromagnetic interactions 
involving odd numbers of sites. Each site may interact with a finite number of other 
sites. As for the case of the Circle Theorem, the dimension and the specific type 
of lattice do not come into play. It is shown that: (i) For systems of the first type, 
there is an interval of the real H axis, (-C(T), C{T)) (with C{T) > that depends 
on the strength of the ferromagnetic couplings) outside which there is no phase 
transition. C{T) goes to zero as T ^ 0, and hence the interval shrinks to H = 0. 
(ii) For systems of the second type, when the interaction involving odd numbers 
of sites is antiferromagnetic, there is an interval on the real H axis (— oo,— C(r)) 
in which no phase transition occurs; when the interaction involving odd numbers 
of sites is ferromagnetic, there is no phase transition for real H G (C(T), oo). Here 
again C(T) is a positive parameter, depending on temperature and on the inter- 
action constants, that goes to zero in the zero-temperature limit. These general 
results put clear constraints on the region of the plane temperature-magnetic field 
where one is entitled to look for a possible phase transition for such systems with 
mul;ti-spin interactions. 
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(B) Ising systems with antiferromagnetic interactions or, equivalently, lat- 
tice gases with repulsive interactions. Various Ising lattice systems with anti- 
ferromagnetic interactions have been considered in Refs. 22 (a Tonks gas of hard rods 
and a gas with a weak long-range repulsion), 31 (Tonks gas), 25 (a lattice gas with a 
hard-core repulsion extended over several lattice sites), 103 (a hard-core lattice gas 
on an M X 00 square lattice), 26, 170 (Ising antiferromagnet), 59 (a two-sublattice 
antiferromagnet), 105 (Ising ferromagnet on a bipartite lattice), 44, 45 (Ising lat- 
tices with various combinations of ferromagnetic and antiferromagnetic couplings), 
46 (Fisher and Temper ley lattice gases), 47 (one-dimensional XY model), 106 (an- 
tiferromagnetic Husimi-Temperley model), 63 (lattice gases with all the zeros of 
the partition function on the real negative semi-axis of the fugacity), 104 (various 
lattice systems of rigid molecules), 107 (various Ising antiferromagnets), 109 (nu- 
merical results on finite square-lattice Ising antiferromagnets) . Different approaches 
are used, that lie from the very abstract rigorous mathematical ones, to approximate 
analytical ones in solving the inverse problem (Sec. 4), and to numerical methods. 
The generic emerging result is that, due to the presence of the antiferromagnetic 
coupling (or, equivalently, of the repulsive part of the lattice gas interparticle in- 
teraction), a part of the Yang-Lee zeros lie on the real negative semi-axis of the 
fugacity z (for finite systems, as well as in the thermodynamic limit, when even- 
tually these zeros occupy densely a whole portion of the negative semi-axis). The 
other details of the distribution of the zeros are, of course, specific to the details 
of the considered model (i.e., range of interactions, whether the coupling is mixing 
ferromagnetic and antiferromagnetic features, etc.). The locus of zeros and even the 
density of zeros were sometimes computed analytically. 

(C) Degenerated Ising spins. In Ref. 110 it is discussed the case of a a one- 
dimensional, nearest-neighbour interactions, spin-1 Ising system that has a spin 
degeneracy, i.e., for which each spin variable S can take either of the values values 
5 = 1,0, —1 with a certain weight (degeneracy) ^(l), 5(0), 5(— 1)- It is shown that 
the corresponding Yang-Lee zeros do not lie, in general, on the unit circle in the 
complex-fugacity plane. The conditions under which one recovers the Circle Theo- 
rem are discussed, see also Refs. 50, 111, 112. 

(D) Van der Waals gas. Several papers were devoted to this prototypical model 
of non-ideal gas, that can be seen as the continuum limit of a lattice gas with weak, 
long-ranged attractive forces, see Ref. 25. References 23, 24, 27, 28, 29, 30 addressed 
the problem of the distribution of Yang-Lee zeros, showing that: (i) For infinite tem- 
perature the zero distribution is located on part of the real negative semi-axis of 
the fugacity; (ii) With decreasing temperature the distribution branches off the real 
axis, circumventing the origin symmetrically, on both sides; (iii) Below the critical 
temperature Tc the distribution forms a closed curve around the origin, with a di- 
ameter decreasing exponentially to zero as T ^ 0; (iv) An additional tail of the 
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distribution remains on the negative real axis, but with a density of zeros going 
hnearly to zero as T — > 0. This tail is connected to the repulsive-core singularity of 

the van der Waals interaction potential, see Ref. 113 for a discussion. The density 
of zeros was computed analytically in some limiting cases, e.g., in the vicinity of the 
branching points. 

(E) ± J Ising model for spin glasses. Reference 114 discusses a numerical tech- 
nique for obtaining the distribution of Yang-Lee zeros of the disordered symmet- 
ric ±J Ising model (i.e., an Ising lattice for which the coupling constant between 

two neigbour sites is chosen at random as J or —J), in two and three dimen- 
sions. This distribution of zeros, which determines the analytical properties of the 
configuration-averaged free energy, is the superposition of the zeros of the partition 
function corresponding to each configuration of the coupling constants. The zeros 
are not distributed on smooth curves, and it seems that in the thermodynamic limit 
they occupy (densely or not, the question is not clarified yet) a whole region of the 
complex- 2; plane. The computation of the corresponding configuration-averaged free 
energy allows to approach some problems of this short-range spin glass model - like 
the existence of the Griffiths singularity (a non-analytic behaviour in the param- 
agnetic phase of a random or diluted Ising ferromagnet, see Ref. 115) and of the 
Almeida-Thouless transition line, Ref. 116. 

9. Yang- Lee zeros and the critical behaviour 

For simplicity, the considerations in this section address the particular case of a 
system for which the Circle Theorem is valid. However, they can be adapted, with 
some modifications, to more general locations of Yang-Lee zeros. For example, in- 
stead of the Yang-Lee edge angle 9o in Sec. 9.2, one should consider, in general, the 
smallest distance between the ensemble C of zeros and the critical point. 

9.1. Density of zeros near criticality 

In the vicinity of a critical point the thermodynamic quantities, as well as the 
correlation functions exhibit power-law behaviours associated with a set of critical 
exponents. For example, for an Ising model one has: 

• For the (zero field) specific heat 



c{t) 



t 



—a 



for t > 



(44) 



and 



c{t) 



\t\ 



for t < ; 



(45) 



• For the spontaneous magnetization 

m{t)^\tf, fort<0; 



(46) 
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• For the isothermal susceptibiUty in zero field 

X{t) ~ t-'', for i > (47) 

and 

X(i)--|tr^', fort<0; (48) 

• For the field dependence of the magnetization at criticality 

m{h,t = 0) ] (49) 

• For the correlation length 

^(i) ~ t-", for i > (50) 

and 

^(t) ~ fort<0; (51) 

• Finaly, for the long-distance spin-spin correlation function at criticality : 

G2ir,t)-^:^, (52) 

where t = (T/Tg— 1) is the reduced temperature (T^ being the critical temperature). 
The nine critical exponents a, a', (3, 7, 7', 6, u, v', and r] are not independent, and 
the scaling laws reduce the number of independent exponents to two."'^"'^^ 

In the framework of the Yang-Lee theory, all the above quantities can be ex- 
pressed in terms of the density of zeros, e.g., \{9, t) for an Ising-like system for which 
the Circle Theorem is valid. In the critical region, \{6, t) should take a particular 
scaling form in order to reproduce the above behaviours. A large body of work has 
been devoted to the study of this scaling form, see Refs. 52, 118, 119, 120, 121, 122, 
123, and the main result is that: 

\{e,t) = \tfij{et-^'^+^'^), (53) 

in the region of the unit circle 6o(t) < \9\ ^ tt where the density of zeros is different 
from zero (see Eq. (38) and above it); i^ix) is a positively-defined function, that is 
finite and nonzero at x = 0^. Thus the critical exponent j3 can be extracted from 
the scaling behavior of A(0, t) in the vicinity of the critical point. 

9.2. Yang-Lee Edge Singularity 

Let us consider again a model for which the Circle Theorem is valid, as for example 
an Ising-like model. For a temperature i > 0, a gap will show up in the density of 
zeros X{9,t). No zeros are present in an interval [—9o{t), 9o{t)] or, correspondingly, 
in a purely imaginary external field interval [—iHo, +iHo] (with Hq = kBT9o/2 in 
the appropriate units, see Sec. 5). This gap reduces to zero when the temperature 
descends to the critical value and below. Kortman and Griffiths "^^ were the first to 
point out the interest to investigate the behavior of X{9, t > 0) for 9 close to the 
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points ±0o{t), called the Yang-Lee edges. A large literature addressed this question, 
see Refs. 113, 124, 125, 127, 128, 129, 130, leading to the following picture. The 

edge singularity should be regarded merely as an usual critical point occuring with 
a purely imaginary field iH. It turns out that 

A(6'. t > 0) - |(9 - 6'o(t)r, for 6 9Q{t)+ . (54) 

Renormalization group arguments demonstrate that the appropriate model describ- 
ing the edge singularity is the so-called (j)^ field theory. The exponent a is simply 
related to the thermodynamic exponent S: 

a = (55) 

The upper critical dimension (above which mean-field predictions are correct) is 

du = 6. Moreover, the above power law behavior is valid for all temperatures T > Tc, 
and also for all 0{n) models. Indeed, the presence of the magnetic field breaks the 
original 0{n) symmetry and thus the behavior becomes Ising-like. Note moreover 
that the critical theory of the Yang-Leo edge in two dimensions corresponds to a 
rather simple realization of conformal symmetry. 

Therefore the behaviour of the density of zeros near the critical edge provides 
a second independent critical exponent S. The knowledge of /3 and S, and of the 
scaling laws characterizes completely the critical behaviour of the system. 

10. Finite-size scaling and Yang-Lee zeros 

As we have seen above, the density of zeros near the critical edge and at criticality 
contains all the information concerning the critical behavior of a given system. 
However, an analytic expression for the density of zeros can be obtained only for 
some particular simple systems. For more complicated cases, the density of zeros has 
to be computed numerically on finite systems. It is thus important to understand 
the role played by the finiteness of the systems. 

The finite-size scaling theory provides a powerful tool in interpreting finite-size 
systems data. The theory was first developped for continuous phase transitions 
based on phenomenological arguments. A more modern presentation is based on 
a renormalisation group approach. The key hypothesis is based on the premise 
that only two different scales matter, namely: ^, the correlation lenght in an infinite 
system, and L, the characteristic linear extent of the finite-size system. The case of 
finite-size effects in first-order phase transition is somehow more subtle as discussed 
by Fisher and Berker, Ref. 134. However, scaling at a first-order transition can 
be treated in a similar manner as the one for a second-order transition with the 
temperature and magnetic anomalous dimensions, yt and yh, assuming the maximal 
values yt = yh = d, where d is the dimension of the system. 

For a finite-size system of size L the singular part of the magnetization has the 
form 



m{t,h,L) = L-'^+y^ m{tLy\hLy>^) 



(56) 
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where h = 2H/kBT. Then the density of zeros scales as 

\{0, t, L) = L-'^+VH x(6Ly^ , tiy' ) . (57) 
At the critical temperature t = one has 

x,{0,L) = L-''+y'^Xc{eLy^), (58) 

which implies for the infinite system 

XciO) r-. \e\'/' . (59) 

For t > 0, the angle of the Yang-Lee edge (see the section above) for finite 

systems, 0o{t,L), scales as 

eo{t,L) = L-y^eo{tLy^), (eo) 

leading, for L ^ oo, to 

00 ~ . (61) 

Finite-size scaling for the density of Yang-Lee zeros A(^^, t, L) has been dis- 
cussed by Alves et al^'^^ for the three-dimensional Ising model, by Kenna and 
Lang 194,201,196 Jqj. ^jjg four- dimensional Ising model and for the 0{N) theory at 
the upper critical dimension '^^^ . Creswick and Kim I35,i97 j^g^yg analyzed the critical 
properties of the two-dimensional Ising model by computing exactly the partition 
function of systems of sizes 4 < L < 10 and by extrapolating the data using a 
Bulirsch-Stoer (BST) algorithm. Janke and Kenna i^^'^o^ have recently devel- 
oped a novel numerical technique which allows to determine both the latent heat 
in the case of a first-order transition and the specific heat exponent a in the case 
of a second-order transition. The main point is the study the so-called cumulative 
distribution of zeros, that is defined as: 

Ai(r,t) = / \{s,t,L)ds. (62) 
Jo 

This approach has been successfully applied to the study of 2d and 2>d Ising models, 
= 2, g = 10 Potts model (see Sec. 13), and lattice gauge theories. 



11. Pirogov-Sinai theory and Yang-Lee zeros for first-order pheise 
transitions 

A recent series of papers, see Refs. 137, 138, 139, 140, 141, addresses the problem of 
a more direct connection of the Yang-Lee zeros with the (complex extension of the) 
free energy and the discontinuities of its derivatives, for general ^ lattice spin models 
in a magnetic field, that present first-order phase transitions. The departure point is 
the Pirogov-Sinai formalism for the phase diagram of the system, see Refs. 142, 143. 



"General" means here systems with no particular symmetry properties like, e.g., the spin-inversion 
symmetry invoked in deriving the Circle Theorem 
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Its basic hypothesis are the existence of a finite number of ground states (which are 
the thermodynamic phases) of the system, and of the availabilty of an appropriate 

contour representation, i.e., the fact that the partition function can be written as 
the sum over the partial partition functions of each phase. Roughly speaking, this 
corresponds to the fact that, for a system with V possible phases {V = 1,2, 3,...), 
the probability distribution for the order parameter (the magnetization) conjugated 
to the external fi(^ld H has V peaks; the partition functions corresponding to the 
different phases can be evaluated from the position and width of these peaks (which 
depend on the control parameters temperature and magnetic field H). The total 
partition function is then expressed as 

V 

ZN{z) = J29if^M-yMz)/kBT] + O(exp(-L/Lo)) , (63) 
1=1 

where N is the number of spins; /; (z) is the complex extension of the thermodynamic 

free energy per unit volume of phase Z, as a function of the complex fugacity z = 
ex.p{—2H / ksT); gi is the degeneracy of phase l;V = L'^ is the volume of the lattice 
of linear extension L and dimension d; and Lq is of the order of the correlation 
length. Let us suppose here, for simplicity, the absence of triple or higher-order 
coexistence points (the theory can be extended to include these cases, too), i.e., for 
any values of the control parameters {T, H) the probability distribution function of 
the order parameter has no more than two peaks. Then by inspection of Eq. (63) 
one realizes that the zeros of the partition function in the complex-z plane arise, up 
to order C'(exp(— L/Lq)); from a destructive interference of pair of terms of the smn, 
qi exp{—Vfi/kBT) and Qp ex.p{—Vfp/kBT), / ^ p = 1, V. So, each zero of Zn{z) 
in the complex- fugacity plane lies in a vicinity of the order 0(exp(— L/io)) of a 
solution of the equations corresponding to this condition of destructive interference: 

Tlefiiz) - {kBT/V)lngi = Tlefp{z) - {kBT/V)\ngp , 
V/kBT[Imfi{z)-Imfp{z)]=TTmod27r, l^p=l,...,V. (64) 

In the thermodynamic limit, the zeros concentrate asymptotically on the phase 
coexistence curves TZefi{z) — TZefp{z), and the corresponding local density of zeros 
is given by p{z) = {2TTkBT)~^ d[fi{z) — fk{z)]/dz. Practically, this means that the 
zeros of the partition function can indeed be expressed in terms of the complex free 
energy and of the discontinuities in its derivatives. References 137, 139 illustrate 
the application of this rather abstract-looking method to three examples, the low 
temperature Ising and Blume-Capel, and the g'-state Potts model in the limit of 
large q. 

12. Fisher zeros of the canonical partition function 

In a review article of 1965, Ref. 4, Fisher suggested the analysis of the phase tra- 
sitions in the frame of the canonical ensemble, through the study of the distribution 
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of the zeros of the canonical partition function in the plane of the complex temper- 
ature - i.e., what we call today Fisher's zeros. Indeed, in view of the equivalence of 

the various statistical ensembles in the thermodynamic limit, one should be able to 
equally characterize a phase transition in any of these ensembles. As in the grand- 
canonical case, the Fisher zeros will lie off the real positive temperature semi-axis 
for a finite system, but may close up on this semi-axis in the thermodynamic limit; 
the point where a limiting line of zeros cuts the semi-axis will locate a critical 
point. Fisher illustrated this idea in Ref. 4 on the two-dimensional square-lattice 
Ising model, with isotropic nearest-neighbour interactions, in the absence of the ex- 
ternal magnetic field, whose partition function is known analytically since Onsagcr, 
see Ref. 144. Fisher found that the zeros of the canonical partition function lie on 
an unit circles in the plane of the complex variable v = smh{2 J /ksT), where J > 
is the ferromagnetic coupling constant, and T is the complex temperature. These 
circles cut the real-u axis in f = ±1; the point v = +1 corresponds to the ferromag- 
netic transition, and v = —1 corresponds to the antiferromagnetic phase transition 
point (in agreement with the ± J symmetry property of the square lattice in zero 
magnetic field). 

The concept of Fisher zeros (both in the absence and in the presence of a mag- 
netic field), which seems to follow closely that of the Yang-Lee zeros studies, was 
rapidly set-up in its general frame, including the relevance of the thermodynamic 
limit, the location of the critical point, Refs. 4, 145, 146, and the characterization 
of the transition, see Refs. 13, 14, 15, 16, 17, 18, 19. The (often tacit) underlying 
assumption of these abstract studies on the characteristics of the transition is that 
in the thermodynamic limit the Fisher zeros fall on smooth, complex-conjugate 
curves, at least in some vicinity of the critical point. In particular, provided that 
this is indeed the case, Itzykson et al, see Ref. 147, proved two very strong general 
statements for Ising systems, namely: (i) These complex-conjugate curves in the 
absence of magnetic field, in the vicinity of the critical point, form an angle with 
the real temperature axis that is an universal known function of the critical expo- 
nent a (of the specific heat) and of the ratio of the specific heat amplitudes below 
and above the critical point, (ii) The angle (with respect to the real temperature 
axis) at which these zeros curves (in the vicinity of the critical point) depart when 
a real, weak magnetic field is switched on is an universal function of the critical ex- 
ponents /3 (of the spontaneous magnetization) and S (of the relation between field 
and magnetization at the critical temperature). 

However, despite some exceptions (like, e.g., the above-cited Onsager-Ising 
model), the above assumption on the location of Fisher zeros is rarely fulfilled, 
see, e.g., Ref. 148 for a discussion of this point. The rule seems to be that the 
Fisher zeros in the thermodynamic limit occupy densely some whole areas in the 
complex-temperature plane, delimited by curves on which the density of zeros may 
diverge. Nothing like a Circle Theorem or an Asano-Ruelle type of argument exists 
for Fisher zeros, because the canonical partition function with respect to a variable 
of the type v = sinh(2 J/fc^T) (where J is a typical interaction constant of the 
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model, and T is the complex temperature) does not have the polynomial structure 
the grand-canonical partition function has with respect to the fugacity z. In fact, 
there are three elements that render the study of Fisher zeros more complicated 
than that of Yang-Lee zeros: 

(i) The first element is that the location of the Fisher zeros is strongly depen- 
dent on the details of the system (dimension, interaction between the components 
of the system, etc.). For example, Fisher zeros for the Onsagcr-Ising model have 
completely different location than the Fisher zeros of the corresponding anisotropic 
ferromagnet. This general feature is, of course, in relation with the fact that the 
critical temperature Tc is a quantity that varies from one system to another. 

(ii) In general, the zeros for a finite system are not necessarily located in the region 
of the complex plane occupied densely by the Fisher zeros in the thermodynamic 
limit. Sec, for example, the case of a mean-field Ising ferromagnet, in Refs. 149, 150, 
151, 152, 153, 101, where the location of the Fisher zeros in the thermodynamic limit 
is interpolated from that corresponding to larger and larger finite systems. 

(iii) The third point is that rather always the location of the Fisher zeros, even in 
the thermodynamic limit, is not well-behaved, i.e., it does not correspond to smooth 
curves in the complex-w plane, but rather to a combination of curves and densely- 
covered regions of the plane (like, for example, for the very simple example of an 
anisotropic ferromagnet on a square lattice, see, e.g., Ref. 148). 

These elements explain the scarcity of exact analytical results on the Fisher ze- 
ros (as compared to the Yang-Lee zeros), the predominance of numerical results, 
and thus the importance of finite-size scaling arguments in the study of the char- 
acteristics of the transition by mean of Fisher zeros. Without entering into details, 
let us give a few examples of the systems that were studied in the literature: 

(A) Ising models with nearest-neighbour ferromagnetic interactions, in- 
cluding: 

(i) Two-dimensional finite m x n isotropic square lattices, in zero magnetic field, 
Ref. 154, 155; see also Ref. 156 for a discussion of the role of boundary conditions. 
Reference 157 generalizes the result of Fisher ^ on the circular location of the zeros 
to m X 00 rectangular lattices with (asymmetric) self-dual boundary conditions, and 
closed form expression of the density of zeros are obtained for m = 1, 2. 

(ii) The Onsager-Ising model in an external, symmetry-breaking magnetic field of 
complex value H/ksT = ±i7r/2, whose exact partition function in the thermody- 
namic limit was first presented by Lee and Yang in Ref. 12, was approached from 
the point of view of the complex-temperature zeros in Refs. 158, 159. The density 
of Fisher zeros is discussed in Ref. 160 (also for the cases of triangular, honeycomb, 
and Kagome lattices). 

(iii) The corresponding problem in a real nonzero magnetic field was approached 
numerically and analytically (with low-temperature series expansions) in Refs. 170, 
159. It was found that the density of Fisher zeros diverges at a non-physical critical 
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point, the Fisher edge singularity that can be considered as the equivalent of the 
Yang-Lee edge singularity (see Sec. 9), as discussed also in Ref. 161. 

(iv) The loci of the Fisher zeros in the thermodynamic limit of two-dimensional 
systems on honeycomb, triangular, diced, and Kagome lattices, in the absence of 
the magnetic field, were described in Ref. 162. Reference 160 discusses the corre- 
sponding density of the Fisher zeros. 

(v) Two-dimensional Ising systems on square, triangular, and honeycomb lattices 
are investigated both in the isotropic and anisotropic cases, for finite systems, and 
in the thermodynamic limit, in zero magnetic field, in Refs. 148, 163, 164, 165, 
166, 167, 168. It is concluded that, in general, the Fisher zeros occupy densely the 
complex-i; plane in the thermodynamic limit, and the corresponding density is also 
discussed. 

(vi) Fisher zeros for finite three-dimensional isotropic, cubic-lattice Ising ferromag- 

nets, in zero magnetic field, as well as estimates (based on finite-size scaling argu- 
ments) of their densities, were discussed, at various levels, in Refs. 169, 170, 171, 
172. 

(B) Higher-spin Ising model. Ising models of spin S=l, 3/2, 2, 5/2 and 3 have 
been calculated in Ref. 173 for a two-dimensional Ising model on a square lattice, 
for various lattice sizes, in zero magnetic field. The temperature-dependence, as well 
as the asymmetry in the specific heat amplitude above and below the critical point 
are estimted from the Fisher zeros distribution. 

(C) Mean-field spin models. The problem of the Fisher zeros for mean-field 
1/2-spin Ising models - including the relation with critical exponents, finite-size 
scaling aspects, the use of a partial distribution function in order to characterize 
the metastable behavior, the relationship with long-range interaction systems - was 
approached in Refs. 43, 149, 150, 151, 152, 153, 101. The exact locus of the Fisher 
zeros for the mean- field spherical model was obtained in Ref. 174. See also Ref. 175 
for an Ising-type mean-field model for a polypeptide with helix-coil transition. 

(D) Antiferromagnetic Ising models. References 170, 176 obtains numerically 
the location of the Fisher zeros for a finite Ising antiferromagnet on a square lattice, 
for different values of the external magnetic field and difii'erent extensions of the lat- 
tice. It is inferred that in the thermodynamic limit these zeros have an accumulation 
point on the real axis, corresponding to the antiferromagnetic critical point. 

(E) Heisenberg ferromagnets. References 177, 178, 179 discuss the thermo- 
dynamic functions (such as energy, specific heat, magnetization, and susceptibility) 
that have been computed numerically, using Fisher zeros distributions, for finite 
three-dimensional Heisenberg models with different anisotropy constant parame- 
ters, including both ferromagnetic and antiferromagnetic couplings. 
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(F) Ising models on hierarchical lattices. The location, fractal nature, and 
density of Fisher zeros for Ising models on various such lattices are studied in 
Rcfs. 180, 181 (diamond hierarchical lattice), 182, 92 (Caylcy trees). 89. 183, 184 
(Sierpinski gaskets and 3- and 4-simplex lattices). We also refer the reader to a 
short, comprehensive review of Itzykson et al. on this subject, Ref. 185. 

(G) Ising models on aperiodic lattices. References 93, 94, 95, 96 discuss the 
Fisher zeros for several such one- and two-dimensional Ising aperiodic lattices, with 
coupling constants generated according some predefined algorithms, see also Sec. 8. 

(H) ± J Ising spin model. In Refs. 186, 187 there are presented numerical results 
on the corresponding Fisher zeros for two- and three-dimensional lattices of differ- 
ent sizes. More extended simulations in Ref. 188 seem to indicate a fractal nature 
of the distribution of zeros. 

(I) Random energy model. The random energy model is one of the simplest, 

exactly-solvable, disordered model that contains some physics of the spin glasses, 
see Ref. 189. It is equivalent to a spin model with multispin interactions exhibiting 
a low-temperature spin glass phase. Its Fisher zeros were studied numerically and 
analytically in Refs. 190, 191, and they were find to occupy densely lines and ex- 
tended areas of the complex-temperature plane. 

Most of the above-cited references invoke finite-size scaling arguments in order to 
infer the properties of the Fisher zeros (location zone, position of the critical point 
of the transition, density, etc.) in the thermodynamic limit from those obtained 
numerically for finite-size systems. A long list of references concentrate precisely on 
the finite-size scaling theory for Fisher zeros (and, sometimes, they address also the 
same problem for the Yang-Lee zeros), see Refs. 192, 193, 194, 196, 197, 198, 199, 
200, 201, 202, 203, 204, 205, 206. 

13. Potts model: what to complexify ? 

The Potts model (so named after R. B. Potts who first studied it in its 1951 Ph.D. 

thesis at Oxford) is a generalization of the Ising model to more-than-two-componcnt 
spins, and, due to its complex behavior, as well as the various experimental realiza- 
tions, has been a subject of intense research during the last thirty years, see Ref. 207 
for a review. The q'-state Potts model consists of g-components spins, i.e., spins that 
can point in the q symmetric directions of a hyperthetraedron in (q—l) dimensions, 
described by the unit vectors s", a = 0, ■■■,q — 1; these spins are placed on a lattice, 
and the corresponding Hamiltonian Ti. of the system in the presence of a magnetic 
field H is expressed as: 

H = -eJ2S{si,Sj)-Hj2S{si,s'^). (65) 

{ij) » 
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Here Sj is the direction of the i-th spin as discussed above; s is the couphng constant 
(which is ferromagnetic if £ > 0, and antiferromagnetic if £ < 0); (ij) indicates 
nearest-neighbour pairs; S is the Kronecker delta; and is a fixed direction amongst 
the q possible ones (one says that the field H is coupled to the spin-state Q). The 
model can be modified to include multi-site interactions. 

Using its representation in terms of a graph with coloured edges, Potts model 
can be generalized to arbitrary, non-integer values of q. It is connected with other 
interesting models of statistical mechanics. Potts model with q = 4 is known as 
Ashkin- Teller model, Ref. 208; the case g = 2 corresponds to the usual Ising model; 
q = 1 model is in correspondence with the percolation process; q = 1/2 can be 
related to a dilute spin glass model; finally, the q = limit relates to the Kirchhoff 
resistor network problem. 

The nature of the order-disorder phase transition in Potts model dc^pends - 
besides on the ferromanetic or antiferromagnetic nature of the model both on the 
dimension d and on the value of q. For example, for the ferromagnetic model in the 
absence of a magnetic field, for d = 2 there is a second-order phase transition when 
g = 2, 3, 4, while for q > 4 the transition is first-order. There exists a dimension- 
dependent critical value of q, qc{d), above which the transition is mean- field like 
- either of the first order if g > max(2, qc), or continuous qc ^ q ^ 2. The few 
exactly- known points are gc(2) = 4 (Ashkin- Teller), gc(4) = 2 (Ising), gc(6) = 1 
(percolation). Besides the q = 2 (Ising) case, the q = 3 and q = 4 models have a 
rather well-studied phase diagram, which proved to be dimension-dependent. 

The location and characteristics of the phase transition(s) - that are known 
analytically only in few cases can be studied, essentially numc;rically, in the Yang- 
Lee formalism. The canonical partition function, in the general case, is a function of 
three control parameters: the magnetic field intensity H, the temperature T, and the 
parameter q of the model. One can therefore consider the complex extension of the 
partition function and the corresponding zeros either with respect to the complex 
magnetic field (Yang-Lee zeros), or the complex temperature (Fisher zeros), or, 
moreover, with respect to the complex q - the so-called Potts zeros. 

13.1. Yang- Lee zeros 

A first category of studies were done on the zeros of the partition function in the 

plane of the complex magnetic field, at fixed temperature and real value of the 
parameter q of the Potts model, for various dimensions d and types of lattices. 

The one-dimensional lattice was studied in the references cited below using a 
transfer-matrix technique. It was shown^^^ that for I < q < 2 the zeros lie inside 
the unit circle in the complex fugacity plane, while for q > 2 they lie outside the 
unit circle for finite temperature; q — 2 corresponds to the Circle Theorem for 
the Ising system. The situation g' < 1 is somehow pathologicaP^'^'^^^, since for all 
the temperatures the zeros lie (in part or completely, depending on T) on the real 
axis, and accumulate around a certain point of the real axis. This would mean 
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that one obtains a second-order phase transition at finite temperature, and in the 
presence of a real magnetic field; however, the eigenvalues of the transfer matrix 
- which give the correlation length - are not real. = 1 is a particular limiting 
case. The discussion of a tricritical point and its connection with a special type of 
Yang-Lee edge singularity is given in Ref. 209 for g = 3 Potts model. Reference 212 
demonstrates a connection between the locus of the Yang-Lee zeros and the Julia 
set of the logistic map associated with the renormalization transformation of the 
lattice (and idea already present in the study of Ising models on hierachical lattices, 
see above). Finally, Ref. 214 presents a recursive method (adapted from the theory 
of dynamical systems) for the computation of the partition function zeros. 

The two- and three-dimensional models on regular lattices are discussed in 
Refs. 216, 217, 211, and they are shown^^^ to have the same type of behavior 
with respect to q as the one-dimensional lattices. Yang-Lee zeros on recursive Bethe 
lattices are studied in Ref. 218 for noninteger values of q. 

13.2. Fisher zeros 

As already mentioned in Sec. 12, the Fisher zeros are much more sensitive than the 
Yang-Lee zeros to the details of the interactions (type of lattice, dimcnsionnality), 
and also to finite-size effects (like boundary conditions), and this explains the large 
number of (essentially numerical) papers devoted to the study of various specific 
models. 

For the one-dimensional lattices with arbitrary q it was shown in Ref. 210 that 
the zeros in the complex field plane can be related, through a kind of duality trans- 
formation, to the zeros in the complex temperature plane. Potts models on square 
lattices are discussed, e.g.. in Refs. 219, 220, 108, in the absence of a magnetic field 
and for various values of q. In the region where TZe[cxp{e / kBT)] < 1, Fisher zeros 
lie on the circle \ exp{£/kBT) — 1| = y^, while outside this region the zeros are 
no longer ocuppying this circle and lead to the existence of (nonphysical) points at 
which the magnetization and susceptibilities diverge - the Fisher edges. The prob- 
lem of the Fisher edges is also discussed in Refs. 221, 222, 161, 223, 214, both in 
the absence and presence of an external magnetic field. 

Finally, a large body of literature addresses these problems on various other 
lattices - strips of square lattices^^'*'^^^, honeycomb, triangular, and Kagome 
lattices^^^'^^'^'^^*, and recursive Bethe lattices.^^^'^^^ 

13.3. Potts zeros 

As already mentioned, Potts model can be extended to arbitrary, continuous values 
of q (through the so-called Kastclcyn-Fortuin^^'' representation). Thus the canon- 
ical partition function becomes a function of the parameter q, and this gives the 
motivation for the analysis of its zeros in the plane of complex q - the Potts zeros. 
The physical relevance of such an approach may be somewhat unclear at first sight. 
However, in Ref. 222 it is shown that one can find a similarity between the scaling 
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property of complex-q zeros and the den Nijs^"^" expression for the thermal critical 
exponent; in Ref. 223 it is shown that a Circle Theorem holds for a certain range of 
intensity of the applied magnetic field, while for other ranges the Potts zeros lie on 
the positive real axis (without necessarily being accumulation, i.e., phase-transition 
points). A Potts edge exists, analogous to the Yang-Lee and Fisher edges. Several 
types of lattices are considered in Refs. 214, 224, 225, 228. Finally, ID long-range 
interaction and the mean-field cases are approached in Ref. 231, with a discussion 
of the phase transition point and of finite-size effects. 

While most of the above-cited references address the case of the ferromagnetic 
Potts model, there is a non-negligible amount of work dedicated to the antiferro- 
magnetic case, see Refs. 108, 224, 228, 223, 222, 213. In particular, the canonical 
partition function in the T ^ limit reduces to chromatic polynomials in q, that 
have a relevance in the theory of graphs, Ref. 232. 

PART II: NONEQUILIBRIUM STEADY-STATE PHASE 
TRANSITIONS 



14. Generalities 

The study of nonequilibrium systems reveals its importance once one realizes that 
equilibrium in nature is merely an exception rather than the rule, and that most of 
the qualitative, structural changes (like, for example, pattern formation) that one 
encounters in various systems take place under nonequilibrium conditions. 

Prom a macroscopic point of view, roughly speaking, a nonequilibrium system 
is the set of fluxes of various characteristic quantities (e.g., number of particles, 
energy, etc.) inside the system and between the system and its surroundings. Let 
us try to give a more precise, quantitative meaning of the notion of nonequilibrium 
system from a microscopic point of view, in the frame of a stochastic description of 
the system in the corresponding configuration space. Of course, as for equilibrium 
systems, this probabilistic description is imposed by the huge number of microscopic 
degrees of freedom, and by the practical necessity of describing the system in terms 
of a few relevant measiirable quantities. 

A configuration iv of the system represents a set of relevant mesoscopic vari- 
ables, corresponding to the degree of coarse-graining of the adopted description. 
The ensemble of the accessible configurations, according to the constraints imposed 
on the system, form the corresponding configuration space. The system is jumping 
between the configurations according to some prescribed stochastic transition rules. 
The state of the system at a time t is given by the probability distribution P{u}, t) 
associated to the set of possible configurations {oj}. Suppose now that the level of 
coarse-graining is such that the evolution of the system is Markovian (this is the 
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only type of systems that we shall be dealing with in the foregoing). The evolution 
of the probability distribution function can be described by a master equation with 
a set of transition rates between the configurations, {>V(a; <^')}-^ 

dP{co, t) 



J2 [ ^{^' uj)P{uj', t) - W{(J uj')P{w, t)] . (66) 



dt 

Let us consider the long-time behavior of the system. The dynamics can be 
such that no stationary state is reached (e.g., the system may exhibit a limit-cycle 
or a chaotic- like behaviour). On the contrary, suppose that the system reaches 
a stationary state, i.e., the probability distribution of the possible configurations 
becomes time-indcpcndcnt, Ps{uj). The stationarity condition in the configuration 
space is that the total flow of probability into configuration u) is balanced by the 
corresponding outgoing flow. 



[ W{uj' uj)Ps{oj') - W{uj uj')P,{uj) = for all u . (67) 



One realizes that an equilibrium state is a very particular case of stationary state. 
Besides the independence on time of the characteristic quantities of the system, 
there is no macroscopic exchange between the system and its surroundings, i.e., no 
flow runs through the system and its borders. Prom a stochastic point of view, this 
corresponds to the detailed-balance condition in the conflguration space, 

W{uj' uj)Ps{oj') = W(a; ^ uj')Ps{i^) for aU uj, J, wi^J , (68) 

i.e., the balance of probability flow between any pair of configurations. 

On the contrary, a nonequilibrium stationary state corresponds to probability 
flow loops in the configuration space, i.e., to the breaking of detailed balance for 
certain configurations w, u', 

W{lo' ^ij)Ps{uj')^W{iJ^uj')Ps{uj). (69) 

These probability flow loops lead to the observed flows of macroscopic quantities 
through the system and its borders. In particular, note that in order to maintain 
a system in a nonequilibrium state continuous exchanges with its surroundings are 
necessary, i.e., the system is necessarily opened. 

One should realize that for such stochastic models the physics is embedded 
in the transition rates, that are chosen on "reasonable" backgrounds according to 
the nonequilibrium phenomena one is wishing to describe. It seems thus that one 
has a high degree of freedom in the choice of the transition rates that lead to a 
nonequilibrium stationary states, although it has been argued recently^^^'^^^'^^^ 
(on the basis of Jaynes' MaxEnt principle extended to nonequilibrium situations) 
that there are rather severe restrictions on this choice. 



^One can consider the more general case when these transition rates are time-dependent but, for 
simplicity, here we shall limit ourselves to the case of constant transition rates. 
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To sum-up, breaking of detailed balance is the general characteristic of stochastic 
systems in a nonequilibrium steady-state. 

15. A partition function for nonequilibrium stationary states 

Phase transitions in nonequilibrium steady states arc somctinics accompanied, as 
in the equilibriTim case, by a breaking of ergodicity of the system, i.c;.. by the fact 
that the asymptotic stationary state of the system is not uniquely-defined, but de- 
pends on its initial configuration. This means that Eq. (67), for certain values of 
the transition rates (depending on the control parameters of the system), admits 
more than one solutions for the stationary probability distribution Ps{(jj), and, cor- 
respondingly, the stationary properties of the system (e.g., the mean values of the 
characteristic parameters and of their stationary flows throughout the system) differ 
from one phase to another. 

The problem we address now is whether there exist any possibility to define a 
partition function for the nonequilibrium systems, that allows for the character- 
ization of a nonequilibrium phase transition between stationary states in a way 
that is analogous to that described in Part I for the equilibrium cases. The answer, 
in a rather general frame, was given in Ref. 10, and we present here the relevant 
arguments, see also Refs. 8, 237, 236, 239. 

For simplicity, define the transition matri:<?^^ [Vl^(a;, w')], whose off-diagonal 
terms are equal to the transition rates W{uj,ix>') = yV{uj' w) and the diago- 
nal terms are W{u),u)) = — J^w' W(w w'), so that the evolution equation (66) for 
the probability density reads 



and it preserves the normalization condition J2oj Pi'^^t) = 1- The stationarity con- 
dition (67) reads 



Suppose now that the control parameters of the system are such that there is a 

unique steady state. Then, according to the Perron-Frobenius theorem, there exists 
a single eigenvalue of the transition matrix that is equal to zero, and the correspond- 
ing eigenvector gives the stationary-state probabilities of the various configurations 
Ps(w). All the other eigenvalues {Xi} of the transition matrix have negative real 
parts, TZeXi > 0. and each of them corresponds to an eigenvector of the transition 
matrix that relaxes exponentially to zero within a time-scale of the order l/|7?,eAi|. 

Because of the zero eigenvalue of the transition matrix, the linear set of equations 
(71) is underdetermined, and therefore one can obtain from it only the relative 
statistical weights fs{u)) of the various configurations. However, one can define the 




(70) 




for all uj . 



(71) 
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such that 

Psico) = ^ fs{u^) . (73) 

This normalization factor is related to the nonzero eigenvalues of the transition 
matrix, more precisely (up to a multiplicative factor),*^ 

Z=l[{-Xi). (74) 

This is the key relationship that allows to propose Z as a candidate for a nonequi- 
lihrium steady- state partition function that can give informations on the appearance 
of a phase transition in the system. 

Indeed, suppose that one varies the control parameters (and thus the transition 
rates, as well as the corresponding eigenvalues of the transition matrix) so that the 
system is driven towards a nonequilibrium phase transition. A new stationary state 
appears in the system, and, correspondingly, at a phase transition one is expecting 
the appearance of a diverging time scale related to this new slow-mode that sets in. 
This means that one of the eigenvalues Aj that is associated to this mode goes to 
zero as one approaches the phase transition point and thus, by virtue of Eq. (74), 
so does the normalization factor Z. 

It seems therefore meaninghil to try to locate the nonequilibrium phase transi- 
tion points through the analysis of the zeros of the normalization factor Z . However, 
the legitimacy of this approach has to be checked on several known exactly-solvable 
models before trying to use it for further predictions. Several classes of nonequilib- 
rium systems were discussed in the litterature using this approach, and we briefly 
review them here, following for the most the original work of Blythe and Evans, 
Ref. 10. 

Recall at this point that Z is defined up to a multiplicative factor,' that depends 
on the method adopted to solve Eqs. (71); this factor (which is commoun to all 



'^In order to obtain this result, note that the characteristic polynomial of the transition matrix 
reads (up to a multiplicative factor) 



det(A/ - W) 



det(W)-hA^/^H + 0(A2 



Given that det(W) = 0, and = X]w ^"^^ 
det(A/ - W) 



allAj AiT^O 



'Note that such a spurious multiplicative factor may appear also in equilibrium systems, e.g., 
through a uniform shift of the energy scale. 
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the stationary states weights fs) is generally a polynomial in the transition rates 
W(w — > w'), and therefore may introduce additional zeros to Z. However, one 
expects that these zeros are not physically relevant, and this can be checked on the 
models that are discussed below. 

16. Driven diffusive systems 

Most of the studies on the applicability of Yang-Lee theory to nonequilibrium phase 
transitions refer to simple models of systems with diffusion and drift, mostly because 
they were intensively studied in the last years and many results are known for their 
stationary states and phase diagrams. 

16.1. Stochastic single- species models 

Such models belong to the class of the so-called ASEP (for "Asymmetric Simple 
Exclusion Process"). These are interacting-particle systems, with a single species of 
particles with no internal degrees of freedom, and with hard-core two-body inter- 
actions with their nearest-neighbour sites on a one-dimensional lattice. This hard- 
core interaction with excluded volume is expressed by the exclusion condition that 
each lattice site may be occupied by at most one particle. Therefore this class of 
models may be described by a set of occupation numbers {ni, n2, n^}, where 
rifc = 0, 1 is the number of particles on site A; of a lattice of L sites. The particles 
may jump at random, according to some prescribed stochastic rules, by one lattice 
site, with a bias in the right/left transition probabilities. The system is maintained 
in a nonequilibrium stationary state, with a constant particle flow along it, by con- 
tinuous injection and extraction of particles at its two borders. The specificity of 
each model is embedded in the prescribed stochastic dynamics, i.e., in the proba- 
bilities for the changes of the configuration of a pair of neighbouring sites k and 
A:-|- 1, and also in the continuous or discrete (parallel or sequential) character of the 
dynamics, see, e.g., Refs. 240, 241 for an overview. 

The totally asymmetric exclusion process (TASEP) with open boundaries is per- 
haps the simplest exactly solvable model with a nontrivial stationary behavior that 
includes both a first-order and continuous phase transitions, and, in view of its 
prototypical character, we will describe it here in some detail, see also Refs. 242, 
10. 

In its continuos-time variant, in an infinitesimal time interval dt, one can have 
one of the following transitions: a particle may hop to the right with probability 
dt (i.e., with transition rate equal to imity). provided the receiving site is empty; 
or a particle is injected at the left border (leftmost lattice site) with probability 
adt (provided the left border is empty); or a particle at the right border is re- 
moved from the system with probability /3 dt (provided the right border contains 
a particle). Obviously, the injection and extraction rates, a and /3, are the control 
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parameters that determine the stationary state of the system, that can be charac- 
terized through the density profile, i.e., the set of mean occupation numbers of the 
lattice site, (rife), k = 1, ...,L. The stationary flow of particles J = J{a,/3) exhibits 
non-analiticities in the thermodynamic limit of an infinite lattice L ^ oo, as one 
varies a and /3, and these non-analiticities correspond to phase transitions in the 
system. The phase diagram consists of three regions: (i) A high-density phase with a 
density profile that decays exponentially towards the right boundary. In this phase 
the current J = /3(1 — /3) is controlled by the rate at which the particles are removed 
from the system, /3 < min(a, 1/2). (ii) A low density phase, with a density profile 
that has an exponential decay towards the left boundary. The current J = a{l~ a) 
is again controlled by the smallest of the two rates, a < min(/3, 1/2). At the first- 
order transition line a = < 1/2, the current exhibits a discontinuity in its first 
derivative. The system presents a shock separating regions of high and low densities, 
and this is an example of phase coexistence at a noncquilibrium first-order phase 
transition, (iii) The third phase, for a, (3 > 1/2, is called the maximum current 
phase, since the current assumes throughout the constant value J = 1/4, which is 
the largest possible current for any combination of a and /3. The density profile de- 
cays as a power-law from both boundaries, and therefore the transitions from either 
phase (i) or phase (ii) are accompanied by diverging lengtliscales. Moreover, J has a 
discontinuity in its second derivative, and this corresponds to a second-order phase 
transition. 

The normalization factor for the TASEP is 

and the current J = Z^-i/Zl. In the thermodynamic limit In Zl ~ —Lin J, i.e., 
In J plays the role of a free energy density for the nonequilibrium partition function 

Zl. From here it is routine to apply the Yang-Lee procedure to Z^. For example, 
keeping /3 as a parameter, one finds that the zeros of Z^ in the plane of the complex 
variable a lie on a smooth curve, and accumulate at the phase transition point Uc 
on the real positive semiaxis when L is increased towards the thermodynamic limit. 
Depending on the value of /3, one has: (i) If /3 < 1/2, the curve of zeros passes 
smoothly, with a nonzero density of zeros, through the transition point ac = (3, 
corresponding to a first-order phase transition, see Sec. 3; (ii) If /3 > 1/2, the 
density of zeros decays to zero as one approaches the transition point oic = 1/2, 
and the zero- lines approach the transition point at an angle 7r/4 to the real axis 
(meeting at a right angle), thus indicating a second-order phase transition. 

Therefore the example of TASEP seems to indicate that it is legitimate to study 
the zeros of the normalization factor Z . in the plane of one complexified control pa- 
rameter and in the thermodynamic limit, in order to get the location of a nonequi- 
librium phase transition point. Moreover, the distribution of these zeros offers in- 
formations on the characteristics of the transition - exactly as in the equilibrium 
case. 
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The same conclusions hold for the time-continuous partially asymmetric exclu- 
sion process (PASEP), for which the particles have a nonzero probability to jump to 
the left, corresponding to a transition rate q. The PASEP has thus a supplementary 
control parameter q, and one can follow the change in the locus of the partition 
function zeros in the complex-a plane (at fixed value of (3) with changing q. Also, 
one could look for the zeros of Z in terms of complex q (at fixed values of a and /3). 
In particular, at g = 1 (when there is no bias between right and left transitions), 
one transits to a regime of zero-current in the thermodynamic limit, and this phase 
transition should manifest as an accumulation of the complex-^ zeros in the vicin- 
ity of 1, a scenario that has still to be verified (see Ref. 268 for some more 
comments on this phase transition). 

Two types of discrete-time variants of TASEP were discussed recently in the 
frame of Yang-Lee formalism for phase transitions, namely: 

(i) The parallel-update TASEP, see Ref. 247. At each time step particles are in- 
troduced with a probability a at the left border (if this first site is empty) and 
are extracted with probability /3 at the right border (if this last site is occupied); 
moreover, the particles in the bulk can jump to the right with probability p (if the 
right space is empty), or remain still with probability 1— p. This update is applied 
to all particles simultaneously. Here again one can compute exactly the normaliza- 
tion factor Z — Z{a,f3,p), and look for its zeros in the complex-a plane at fixed f3 
and p. The phase diagram is similar to that of the continuous-time TASEP; there 
are three phases, separated, respectively, by two continuous-transition lines at a or 
/3 = 1 — \/l — p, and a first-order transition line at a = /3 < 1 — \/l — p. Thus the 
zeros of Z in the complex-a plane also behave similarly. One can also look at zeros 
of Z in the complex-p plane (at fixed values of a and /3) ; their locus is a cardioid, 
with the transition point appearing as a cusp on the real positive semiaxis. 

(ii) The sublattice-parallel update TASEP with the supplementary constraint of fixed 
average density of particles in the system, p = M/ L (where M is the total number 
of particles on the lattice of length L) was studied in Ref. 246. In this model, at the 
first time step, besides the injection and extraction of the particles at the borders 
with probabilities a and P, respectively, all the particles on odd sites jump to the 
right with probability one (if the site to the right is empty); at the next step, all 
the particles on even sites jump to the right with probability one (again, if the re- 
ceiving site is empty). Then one repeats this two-stage process. The phase diagram 
of this model in the thermodynamic limit contains three phases, namely a high- 
density phase and a low-density phase (as for the continuous-time TASEP), and 
a shock phase, characterized by existence of shocks in the particle-density profile. 
This shock phase has second-order transition borders with the high and low density 
phases, and this is confirmed by the behavior of the zeros of Z in the complex-a 



JThis is a special case of the famous Nagel-Schreckenberg model for traffic flow [K. Nagel and M. 
Schreckciibcrg, J. Phys. I. France 2, 2221 (1992)], in which the TASEP phase transitions can be 
considered as jamming transitions. 



2, 2008 17:43 WSPC/INSTRUCTION FILE Yang-Lee 



YANG-LEE FORMALISM FOR PHASE TRANSITIONS 41 

plane, for various fixed values of /3 and p. 

The problem of finite-size scaling and universality for several variants of discrete- 
time TASEP was addressed recently in Ref. 248. In the thermodynamic limit, 
in general (i.e., whatever the adopted update rule), both first and second order 
boundary-induced transitions between the steady-states are encountered, and one 
derives finite-size scaling expressions for the current J, for the mean local particle 
density, as well as for the zeros of the normalization factor Z in the complex-a plane. 
In particular, simulations suggest that the smallest distance between the locus of 
these zeros and the critical point Uc has a power-law behavior L~^/" , where u (as 
in the equilibrium systems) is the critical exponent of the bulk correlation length 
for the nonequilibrium phase transition. 

16.2. Two-species models 

Actually, there are two slightly different variants of the same model that were stud- 
ied, from the point of view of the Yang-Lee approach, in Refs. 249, 250. As described 
in Ref. 249, the model consists of two types of particles, A and B, of fixed numbers 
Na and TVs, that diffuse in opposite directions on a one-dimensional L-site lattice 
with periodic boundary conditions. Specifically, particles of type A jump to the right 
and particles of type B jump to the left, both with unit transition rates, provided 
that the destination site is empty. Also, if on neighbour sites, the A and B particles 
exchange places, A-\-B — > B-\-A with transition rate q, and B-\-A — > A-\-B with unit 
transition rate. The bias parameter q is thus the control parameter of the system 
(at fixed particle densities pa.b — Na,b/L), and as long as g 7^ 1 the system can 
reach a nonequilibrium stationary state with a stationary nonzero flow of particles 
throughout the system. Arndt^^® studied the phase diagram in the particular case 
of equal particle densities pA = Pb ^ P = M / L in terms of the Yang-Lee formalism. 
Since the normalization factor Z ^iq, p) was not know analytically at that moment, 
Arndt went to the numerical study of a grand-canonical partition function 

EL{z,q)= Yl z''ZL{q,M/L) (76) 

M=0 

in the plane of the complex fugacity z. From these studies of the zeros on finite-size 
systems (with L up to 100), he inferred the existence, in the thermodynamic limit, 
of a first-order phase transition at some point qc{p) > 1, between a mixed phase 
(corresponding to a condensate made of both A and B particles; a block of empty 
sites, with a few residual A ov B particles, occupies the rest of the system) and 
a disordered phase (with uniform density profiles of both species and no spatial 
condensation) . 

However, it was shown later, see Refs. 10, 251, 252, through an exact calcula- 
tion of the grand-canonical partition function in the thermodynamic limit, that no 

phase transition stricto scnsu appears in the system. When one varies the control 
parameters, there appears, instead, a very abrupt increase in the correlation length 
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to an anomalously large value - of the order 0(10''°) sites, which, however, if finite. 
Although very spectacular, this phenomenon does not seem to be extremely rare, 
since it can be argued, sec, e.g., Rcf. 252, that it is generated when the dynamics 
of the model is such that small domains tend to be suppressed in the steady-state 
distribution. 

Thus the extrapolation done by Arndt from numerical results on finite-size sys- 
tems towards the thermodynamic limit is unreliable. One should therefore be aware 
of the possible appearance of such extremely large correlation lengths in the nonequi- 
librium systems, correlations that may lead to a false impression of the existence 
of a phase transition when investigating (numerically) systems with sizes less or 
comparable to the correlation length. 

Note also^° that investigating the zeros of this grand-canonical partition function 
in the complex fugacity plane amounts, from a physical point of view, at placing 
the lattice in equilibrium with a reservoir of particles, and thus (despite the claim 
in the title of Ref. 249) at studying an equilibrium system. 

Reference 250 reconsidered the same model in the exactly-solvable case Nb = 1, 
i.e., of a single impurity present in the system. Using a matrix-product formalism, 
one can compute the corresponding normalization factr ZL{q, p) at fixed density p of 
the A particles. Studying the behaviour of ZL{q, p) in the thermodynamic limit, one 
deduces the existence of two stationary phases of the system, namely, (i) a jammed 
phase for q < 2p, when the impurity provokes a macroscopic shock (error-function 
like) in the density profile oi A particles; (ii) and a power-law phase for q > 2p, when 
the impurity has a short-range effect on the system, such that the density profile 
of As has an exponential behaviour, with a correlation length ~ \\n{qc/q)\~^ that 
diverges when q approaches the critical value q,. = 2p. This result is confirmed by 
the study of the location and distribution of zeros of Zl (q, p) in the complex-g plane 
(at fixed value of p), that predict correctly the second-order phase transition point 
between the two phases, as well as its characteristics. 

17. Simple reaction-diffusion models 

These stochastic models are characterized by the fact that, contrary to the drive- 
diffusive systems discussed above, the number of particles is not conserved by the 
dynamics, i.e., there are some reactions inside the system. One such model, which re- 
sults through a simple modification of the continuous-time ASEP, is the coagulation- 
decoagulation model with open boundaries, see Ref. 243. It consists of particles that 
diffuse, coagulate (two neighbour particles merge into a single one), and decoag- 
ulate (one particle splits into two neighbour particles) preferentially, e.g., in the 
leftward direction. More precisely, in the simplified model in Rcf. 244, they (i) dif- 
fuse to the left and right with rates q and q^^ respectively; (ii) coagulate at the left 
A + A^ A + with rate q and at the right A + A^O + A with rate q~^; (iii) deco- 
agulate to the left + A^ A + A with rate jq and at the right A + ^ A + A with 
rate ^q~^- Here g > 1 to assure the leftward preference of the processes. Besides 
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this, there is an injection and an extraction of particles solely at the left border, with 
rates a and /?. This model, for which the normalization factor Z can be computed 
exactly when a = {q~^ ^ 9 + /?)7. exhibits a first-order phase transition between 
a low-density and a high-density phase. The roots of Z in the complex-q plane 
are studied (for fixed and 7), and Yang-Lee theory is shown to hold, describing 
correctly both the location and the nature of the phase transition. 

A variant of this model with reflecting boundary conditions, instead of the open 
boundaries, was studied in Ref. 244. The model in this case has two control pa- 
rameters, q and 7, at a fixed particle-density p. The system has two second-order 
phase transition points, at qc = l/\/l — p and q'^ = 1/qc = \f\-~p. The reason 
of the existence of these two transition point is simply the invariance of Z under 
the transformation q — > . Note that 7 plays no role in the location of the phase 
transition point. Again the Yang-Lee formalism predicts correctly these results. 

Reference 245 offered recently a unified frame, based on the matrix-product for- 
malism, for the study of the nonequilibrium steady-state phase transitions of three 
families of one-dimensional nonequilibrium models with opened boundaries. Besides 
the TASEP and a generalization of the above-described coagulation-decoagulation 
model, it also discusses the asymmetric Kawasaki-Glauber process (AKGP). The 
latter corresponds to a (totally) asymmetric diffusion of the particles (Kawasaki 
spin-exchange dynamics), combined to death to the left and right A + + 0, 
0-\-A 0-hO, and to branching to the left and right A-\-0 A-\-A, 0-\-A A-\-A. 
These two last processes represent a variant of the spin-flip Glauber dynamics. More- 
over, a steady-state ciirrcnt is maintained through the injection and extraction of 
particles at the borders. For special values of the corresponding transition rates, one 
can estimate the normalization factor Z, and study its Yang-Lee zeros (however, the 
paper does not concentrate further on this neither for the AKGP, nor for TASEP, 
but only for the extended coagulation-decoagulation model). 

18. Directed percolation 

An important type of nonequilibrium phase transitions is represented by the 
absorbing-state phase transitions, which take place when a system, during its evo- 
lution, reaches a configuration (a so-called absorbing state) in which it remains 
trapped for ever. Phenomena of this type are present in a wide variety of physi- 
cal, chemical, biological systems, see Refs. 253, 255, 254 for recent reviews. It has 
been conjectured recently that a large group of such models fall into the directed- 
percolation (DP) universality class (see Sec. 20 below for some comments on univer- 
sality classes in nonequilibrium systems, and Refs. 256, 257 for this conjecture). A 
simple example of such a model is a reaction-diffusion system on a one-dimensional 
lattice, with a competition between the reproduction of particles, i.e., a (symmet- 
ric) decoagulation, A + O^A + AorO + A-^A + A, and the death of particles, 
i.e., a spontaneous decay A 0. This competition is controlled by some external 
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parameter p; if its value is below a certain threshold p < pc, the system reaches an 
absorbing state (empty lattice) with certainty; otherwise the system remains active 
forever in the thermodynamic limit of an infinite-lattice. The corresponding phase 
transition is a continuous one. 

In order to study the applicability of the Yang-Lee theory to this class of nonequi- 
librium phase transition, it is natural to turn to the original directed percolation 
model, see Refs. 259, 260, 261, 10. In order to fix the ideas, let us consider the two- 
dimensional case - although the process can be defined on general d-dimensional 
lattices. The system consists of a finite square lattice with L rows, that has one 
site, the apex O, on the first row, two sites on the second row, etc., up to a total of 
L{L + l)/2 sites. Bonds between the sites are opened with probability j3, and closed 
with probability 1 — p. At the initial moment, a particle is placed on the apex; at 
each time step this particle jumps on a site of the next lower level, provided the 
corresponding bond is opened. The sites connected by a continuous percolation path 
(i.e., by a succession of opened bonds) are said to belong to the same percolation 
cluster. One is interested by the percolation (or survival) probability Pl{p) that at 
least one site on the level L is connected to the apex through a percolation path; 
Pl{p) is obtained by considering all the possible bond configurations. The order 
parameter, defined for an infinite lattice L — > oo, is P^{p), which is the probability 
of having an infinite cluster (i.e., a cluster that goes from the apex all the way 
down the lattice) for a given bond probability p. 11 p < pc, then Poo{p) is zero, i.e., 
there are no infinite percolation clusters in the system (that means that one can 
always, i.e., with probability one, find a row, suHiciently far from the origin, that 
is not connected to the apex through any percolation path. For p > pc, however, 
P(x {p) becomes nonzero, and in the vicinity of the critical point it has a power-like 
behavior, 

Poo{p)-^{p-Pcf . (77) 

with an universal exponent (3. In addition, the DP process is characterized by a 
transversal correlation length perpendicular to the percolation flow, and a lon- 
gitudinal correlation lenght ^|| along the percolation flow, that diverge both when 
p approaches its critical value Pc, 

a-lp-Pcl""" ^II-Ip-PcI-"". (78) 

A continuous transition takes place at pc- Although DP can be defined and simulated 
easily, no analytical solution exists for it (neither for any of the models pertaining to 
the DP universality class) in dimension d lower than the critical dimension dc = 6. 
Therefore, one has to use numerical estimates of the critical exponents. Currently, 
the best estimates^-'^s in = 2 are /3 = 0.27649, iy± = 1.096854, and = 1.733847; 
the percolation threshold depends both on the spatial dimension and the type of 
array, and for d = 2 on a square lattice the best estimation of it is Pc = 0.6447001. 

Attempts to describe this phase transition in the Yang-Lee formalism were done 
in Refs. 259, 260, 261, by considering the zeros of Pl{p) in the plane of complex 
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p, for dimensions L < 15, and looking for their accumulation points when one 
increases L, i.e., when one goes towards the thermodynamic limit. These zeros do 
not lie on a single smooth curve (as it was the case for the other nonequilibrium 
systems discussed above) , but on a sequence of curves that tend to intersect at the 
critical point. Also, the results in Refs. 259, 260, 261 suggest a fractal structure of 
the distribution of these zeros (but this aspect deserves further investigation). Prom 
these numerical results, the authors infer the values of Pc, and (using finite-size 
scaling arguments) also the longitudinal exponent | . 

It should be noted, however, that the partition function Pl{p) has some prop- 
erties that render it qualitatively different from the other nonequilibrium partition 
functions considered above. A first remark 259,262 jg ^^^^^ ^j^g survival probability 
Pl{p) can be represented formally as the partition function of a system of Ising 
spins, with two-spin and three-spin plaquette interactions, on a pyramid-shaped 
lattice. However, some of these three-spin interactions have an infinite energy, and 
thus the usefulness of such a mapping is questionable. Also, one encounters some 
problems with the definition of a thermodynamic limit and of an extensive free en- 
ergy of the system, since, on one hand, Pl{p) does not grow exponentially with the 
size of the system, and on the other hand hmi^oo Pl{p) = for ^ p < Pc, and it 
is not clear how to define a free energy density in this case. This illustrates some of 
the difficulties related to the extension of equilibrium concepts to nonequilibrium 
situations. 

19. Self-organized criticality 

Two recent studies, Refs. 263, 264, addressed the problem of the applicability of the 
Yang-Lee formalism to self-organized criticality (SOC). Roughly speaking, a sys- 
tem that exhibits SOC is a nonequilibrium system (subject to external stationary 
nonequilibrium constraints) characterized by the fact that, in the thermodynamic 
limit, it reaches a stationary state with scale invariance and power-law statistics. 
One should underline that this nonequilibrium state is reminiscent of an equilib- 
rium critical states; however, contrary to these ones, this stationary critical state 
is reached solely through the internal dynamics of the system in response to the 
external constraints, without the fine tuning of the control parameters that is re- 
quired for the equilibrium criticality. There is a huge body of literature devoted to 
this topic, see Ref. 265 for a pedagogical overview of the field. 

The papers of Cessac et al. rise then the legitimate question whether this spon- 
taneous criticality manifests itself in the properties of the distributions of complex 
zeros of a properly defined partition function for the stationary state, i.e., if these 
zeros have an accumulation point when the size of the system is increased to infinity, 
without any manipulation of the parameters of the system. 

The dynamics of SOC systems occurs in avalanche-like events; one can define a 
set of observables (size, duration, etc.) characterizing these avalanches. In a finite 
system of size L such an observable V can take a finite set of values, up to a 
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maximum finite value Nl, V = 0, 1,...,Nl, according to a stationary probability 
distribution Pl{V = n). Let us consider the generating function of this distribution 
(the grand-canonical partition function), 



that is a polynomial of degree Nl in the fugacity z. If the system is SOC, the zeros 
of this polynomial in the complex-^ plane should have an accumulation point on 
the real positive semiaxis in the thermodynamic limit i — > oo for any parameters 
of the system. This reasoning was verified in the above-cited references for various 
finite-size scaling forms of Phin) encountered in the literature on SOC: the zeros 
pinch the real axis at z = 1 as the system size goes to infinity. A scaling theory for 
the Yang-Lee zeros is proposed in this setting, that shows, under specific conditions, 
a violation of the scaling usually observed in equilibrium critical phenomena. 

Note however that the remark we made in Sec. 16.2 with respect to Ref. 249 
applies here, too, namely that considering the generating function amounts 
actually at analysing an equilibrium problem in a grand-canonical statistical en- 
semble. 

20. Connection with equilibrium systems with long-range 
interactions 

The examples in the above sections show a surprising analogy between the equilib- 
rium and the noncqilibrium steady-state phase transitions, as far as the behaviour 
of the zeros of the respective partition function in the plane of the complex con- 
trol parameter is concerned. The question that arises is then whether or not this 
analogy can be pushed further, e.g., if one can find a kind of (general rule for a) 
correspondence between a nonequilibrium steady-state system and an equilibrium 
one. 

As discussed in Sec. 14, a nonequilibrium system is characterized by a breaking 
of detailed balance in the configuration space, which corresponds to the appearance 
of macroscopic fiows throughout the system. This generates effective long-range in- 
teractions and correlations in the system, as well as the emergence of an efii'ective 
criticality in nonequilibrium steady-states, and, in this logic, of universal distribu- 
tion functions for macroscopic quantities, of the corresponding classes of univer- 
sality, critical exponents, etc., see Refs. 266, 267, 8, 255 for a discussion of these 
fundamental issues. 

Therefore, if one is searching for a kind of correspondence between nonequilib- 
rium steady-states and equilibrium systems, the latter ones have to have long-range 
interaction Hamiltonians. The history of searching for these effective Hamiltonians 
(and related problems like, e.g., definitions of an eff'ective nonequilibrium temper- 
ature) is long and tortuous and will not be discussed here. In particular, trying to 
write down an explicit expression of this Hamiltonian does not seem to be a very 




(79) 
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fruitful approach. 

Here we will briefly present a few recent results in this direction, see Refs. 268, 
269, 247, 270. These studies refer to ASEP-like systems, but in our opinion, in view 
of their generality, these ideas carry an important potential for further applications. 
The adopted strategy relies on two key-elements as discussed below. 

(A) The formal definition of "particle numbers" and associated "fugaci- 
ties". Let us consider a Markovian system as described in Sections 14, 15, and the 
corresponding steady-state normalization factor Z„({W(a; — > w')}) (here n desig- 
nates the number of configurations lu accessible to the system; it is also explicitly 
indicated that the normalization factor is a function of the transition rates) . A first 
point discussed in Ref. 268 is that the normalization factor is a polynomial in the 
transition rates W(w co'), with positive coefficients, of degree n— 1, i.e., it has the 
form of a generating function. By the Cauchy- Schwartz inequality it follows that 
its negative logarithm (that we would like to assimilate with a free energy of the 
system), 

F„ = -ln(Z„), (80) 

is a convex function in all its arguments W{uj uj'). One is then tempted, by 
analogy with equilibrium situations, to identiiy formally the transition rates with 
"fugacities" , and to consider the corresponding "particle numbers" 

= -''I aimZrhTl^ (81) 

These are well-behaved thermodynamic quantities, in the sense that they are pos- 
itive and increasing functions of the fugacities for any size n of the system. They 
are, however, linearly-dependent, and some of them may even coincide (because 
some of th(! transition rates may be zero or may be equal between them). In the 
thermodynamic limit of large n 

iVc.,^.' = V{n)p^^^, , (82) 

where V{n) is the "volume" (defined by the leading asymptotic behaviour of the 
normalization) and Pa;,^' are the "densities" . 

One has now all the formal ingredients to discuss a nonequilibrium steady-state 
phase transition by analogy with the equilibrium situations. For example, a first- 
order phase transition would correspond to a discontinuity of the fugacity as a 
function of the density. It is not clear however, in general, that all the known char- 
acteristics of the phase transitions (e.g., diverging correlation length for continuous 
one) are recovered within this formalism. 

There is one more very important point to be made, and which represents a ma- 
jor diff'erence between nonequilibrium systems and equilibrium systems with short- 
range interactions. It is the fact that the "particle numbers" defined above are not 
necessarily extensive quantities. This implies that in the space of the parameters of 
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the system (which are the transition rates), besides the regions where the "particle 
densities" are well-behaved and finite, there may appear regions where Pa>,a>' may 
diverge (and one has to change then the definition of V{n)). The frontiers between 
these regions correspond then to phase transitions that do not have a correspon- 
dent in the short-ranged interactions equilibrium systems. Note, however, that phase 
transitions of this type are know in equilibrium systems with long-range (or nonlo- 
cal) interactions^^''^ , a remark that goes in the sense of the general statements made 
in the beginning of this paragraph. 

(B) The correspondence with equilibrium systems. Once this formal frame 
is set-up, comes the most delicate part, which is to find a way (possibly a systematic 
one) to assign a physical meaning to the normalization factor, "particle densities" , 
etc. in a properly-defined equilibrium statistical problem. The general solution pro- 
posed in Ref. 268, and illustrated on concrete examples in Rcfs. 268, 269, 247, 
270 comes from the combinatorial graph theory, and consists in relating explicitely 
the normalization of a stationary state to the combinatorial problem of counting 
weighted spanning trees on graphs; the weights of the configurations depend on pa- 
rameters which correspond to the transition rates of the nonequilibrium stochastic 
process. This implies directly an interpretation of the normalization factor as a sta- 
tistical mechanics partition sum. Of course, establishing the right correspondence 
is a difficult task, and no general recipe can be given for it. 

A first heuristic application of this approach is given in Ref. 271, where a 
one-dimensional, stochastic, adsorption-desorption nonequilibrium model for in- 
terface growth (the raise-and-peel model) was put in correspondence with a two- 
dimensional ice model with domain- wall boundary conditions (an equilibrium prob- 
lem with nonlocal interactions). Brak and Essam^^^ showed that the matrix rep- 
resentation of the stationary-state algebra of TASEP can be interpreted combina- 
torially as various weighted lattice paths. The normalization factor of TASEP is 
identical to the equilibrium configuration sum of a polymer chain having a two- 
parameter interaction with a surface (the so-called one- transit walk). As shown in 
Ref. 268, the TASEP current and density are simply related to the equilibrium 
densities of the latter model; in particular, the second-order phase transitions of 
the TASEP (see Sec. 16.1) are related to a special type of surface phase transition 
encountered in polymer physics. Later on Blythe et al found that the normal- 
ization factor of the discrete-time parallel-update TASEP can be expressed as one 
of several equivalent variants of two-dimensional lattice constrained-path problems. 
Finally, in Ref. 270 it was shown that the grand-canonical normalization factor of 
the TASEP (that corresponds to the immersion of the lattice in a reservoir of par- 
ticles, with fixed fugacity z) can be computed in a closed analytical form, and thus 
it allows a direct derivation of the asymptotics of the canonical normalization for 
the various phases, and the immediate correspondence with the one-transit walk 
mentioned above. 
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21. Conclusions and perspectives 

In order to explain the origin of singularities of the thermodynamic potentials in 
some first-order equilibrium phase transitions, Lee and Yang proposed to examine 
the behaviour of the zeros of the partition function. Such an idea turned out to 
have far reaching consequences. Indeed, subsequent generalizations of this approach 
showed that from the analysis of zeros of the partition function one can extract a 
lot information about phase transitions in many different systems. When combined 
with finite-size scaling arguments, this method may be used to calculate the location 
of the transition point and even critical exponents. An important point is the fact 
that in this approach we directly refer to the properties of the partition function, 
which is the most basic quantity of equilibrium statistical mechanics. Thus, one can 
apply mathematically rigorous techniques to obtain a number of interesting results 
concerning, e.g., the very existence and the location of phase transitions in various 
systems. 

Recent largo interest in the Yang-Lcc approach is rclatcxi. however, to its further 
extension to nonequilibrium phase transitions. Several examples show that accom- 
plishing such an ambitious task may be feasible in some nonequilibrium systems. But 
the entire approach is much more problematic than in equilibrium situations. It is 
not entirely clear what is the analogue of the partition function in this case. For the 
directed percolation, coalescence of zeros was observed for the percolation probabil- 
ity - that is actually the order parameter of the system. For self-organized criticality 
this effect appears for some generating functions. Perhaps the most promising are 
the results obtained for TASEP models. In this case the normalization factor of a 
steady-state probability distribution seems to play the role of the partition function. 
In some cases one can even establish much closer relations with equilibrium systems. 
But it is not known yet whether such an approach can be applied to other classes 
of nonequilibrium systems. For example, it would be interesting to test this method 
on systems for which the transition rates are not constant, but state-dependent, 
like, e.g., the zero-range processes (ZRP). 

It seems therefore that, after more than fifty years of generalizations and exten- 
sions, the full potential of the Yang-Lee approach is yet to be uncovered. Important 
developments in its application to the study of nonequilibrium phase transitions are 
in progress. 
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